C SOURCES:
C    Simone Paziani, Saverio Moroni, Paola Gori-Giorgi, and Giovanni B. Bachelet.
C    Local-spin-densityfunctional for multideterminant density functional theory.
C    Physical Review B, 73(15), apr 2006.



C*****************************************************************************
      pure subroutine ESRC_PW92_ERF(rho_c, mu, E)
C*****************************************************************************
C   Implemented by E.R. Kjellgren.
C
C   Subroutine generated using Sympy 1.3
C   Generated: March 25, 2019
C*****************************************************************************
      implicit none
      real*8, intent(in) :: rho_c, mu
      real*8, intent(out) :: E
      real*8 :: x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x
     &13, x14, x15, x16, x17, x18, x19, x20, x21, x22, x23, x24, x25, x2
     &6, x27, x28
      E = 0.0d0
      x0 = 9.86960440108936d0
      x1 = 2.0d0/x0
      x2 = 0.693147180559945d0
      x3 = -x2 + 1.0d0
      x4 = 0.682784063255296d0*rho_c**(-0.333333333333333d0)
      x5 = 2.14502939711103d0
      x6 = rho_c**0.666666666666667d0
      x7 = 1/(x5*x6)
      x8 = 0.564189583547756d0*rho_c**(-0.5d0)
      x9 = 0.826307487110758d0*rho_c**(-0.166666666666667d0)
      x10 = x3*(0.194159335344114d0*x4 + 1.0d0)*log(0.5d0*x0/(x3*( 3.259
     &5509194222925d0*x4 + 0.40691300451752932d0*x7 + 1.4187228164796675
     &d0*x8 + 7.2401019343168302d0*x9)) + 1.0d0)
      x11 = mu**2
      x12 = 1.41421356237310d0
      x13 = 4.60115111447049d0
      x14 = rho_c**(-1.33333333333333d0)
      x15 = 3.14159265358979d0
      x16 = rho_c**1.0d0
      x17 = 1d0/x16
      x18 = x17/x15
      x19 = (-0.009578475d0*x18 + 0.0188171922990732d0*x4 + 0.0676322201
     &645715d0*x7 + 1.0d0 + 0.00126674656487366d0*x14/x13)*exp( -0.68361
     &07611867116d0*x4)
      x20 = rho_c**(-1.66666666666667d0)
      x21 = 0.0837224526465878d0*x20
      x22 = x11*x4
      x23 = 5.57236303607474d0*mu*x9 + 1.0d0
      x24 = mu**3
      x25 = 0.148394183600699d0
      x26 = x19 - 1.0d0
      x27 = 2.28942848510666d0*x13
      x28 = 1.63540853354893d0*x27*(-0.0259335011650457d0*x4 + 1.0d0)/( 
     &0.494402081358784d0*x4 + 0.0524148278841779d0*x7 + 1.0d0)
      E = -rho_c*(x1*x10 - (0.000640558526213818d0*mu**8*rho_c**( -2.666
     &66666666667d0)*x10 + mu**6*(0.0108059179628909d0*rho_c**( -2.0d0)*
     &x10 + 0.0167302167879722d0*x20*x25*x26) + 0.0223069557172962d0*mu*
     &*5*x12*x19*x21 - mu**4*( 0.00139418473233101d0*rho_c**(-1.0d0)*x25
     &*(1.63540853354893d0*x27 - x28 - 10.9027235569928d0*x5*(-0.3525213
     &95009435d0*x4 + 0.558025705063192d0*x7)*exp(-0.49698248213959023d0
     &*x4)) - 0.0683590130188305d0*x10*x14 - 0.13157433081915d0*x18*x26)
     & - x1*( x2 - 1.0d0)*log((6.7683429898050367d0*x22 + x23 + 3.392602
     &5578013106d0*x24*x8)/(3.1331792677937806d0*x22 + x23)) + x12*x24*(
     &0.0315054072231411d0*x17*x19 + 0.00111534778586481d0*x21 *(12.0d0*
     &x15*x16*(-4.49737346725955d0*x4 + 0.825481812223657d0*x7 )*exp(-0.
     &28165369188898165d0*x4) + x28*x6)))/(0.508616435555896d0 *x11*x7 +
     & 1.0d0)**4)
      end subroutine


C*****************************************************************************
      pure subroutine D1ESRC_PW92_ERF(rho_c, mu, E, d1E)
C*****************************************************************************
C   Implemented by E.R. Kjellgren.
C
C   Subroutine generated using Sympy 1.3
C   Generated: March 25, 2019
C*****************************************************************************
      implicit none
      real*8, intent(in) :: rho_c, mu
      real*8, intent(out) :: E, d1E(9)
      real*8 :: x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x
     &13, x14, x15, x16, x17, x18, x19, x20, x21, x22, x23, x24, x25, x2
     &6, x27, x28, x29, x30, x31, x32, x33, x34, x35, x36, x37, x38, x39
     &, x40, x41, x42, x43, x44, x45, x46, x47, x48, x49, x50, x51, x52,
     & x53, x54, x55, x56, x57, x58, x59, x60, x61, x62, x63, x64, x65, 
     &x66, x67, x68, x69, x70, x71, x72, x73, x74, x75, x76, x77, x78, x
     &79, x80, x81, x82, x83, x84, x85, x86, x87, x88, x89, x90, x91, x9
     &2, x93, x94, x95, x96, x97, x98, x99, x100, x101, x102, x103, x104
     &, x105, x106, x107, x108
      E = 0.0d0
      d1E(:) = 0.0d0
      x0 = 9.86960440108936d0
      x1 = 2.0d0/x0
      x2 = 1.46459188756152d0
      x3 = 1d0/x2
      x4 = rho_c**(-0.333333333333333d0)*x3
      x5 = 0.194159335344114d0*x4 + 1.0d0
      x6 = 0.693147180559945d0
      x7 = -x6 + 1.0d0
      x8 = 2.14502939711103d0
      x9 = 1d0/x8
      x10 = rho_c**0.666666666666667d0
      x11 = 1d0/x10
      x12 = x11*x9
      x13 = 0.564189583547756d0
      x14 = rho_c**(-0.5d0)*x13
      x15 = 0.826307487110758d0
      x16 = rho_c**(-0.166666666666667d0)*x15
      x17 = 0.406913004517529d0*x12 + 1.41872281647967d0*x14 + 7.2401019
     &3431683d0* x16 + 3.25955091942229d0*x4
      x18 = x0/(x17*x7)
      x19 = x7*log(0.5d0*x18 + 1.0d0)
      x20 = x19*x5
      x21 = x1*x20
      x22 = mu**2
      x23 = 0.508616435555896d0*x12*x22 + 1.0d0
      x24 = x23**(-4)
      x25 = mu**5
      x26 = 0.0837224526465878d0
      x27 = 1.41421356237310d0
      x28 = exp(-0.6836107611867116d0*x4)
      x29 = 4.60115111447049d0
      x30 = 1d0/x29
      x31 = rho_c**(-1.33333333333333d0)
      x32 = rho_c**1.0d0
      x33 = 1d0/x32
      x34 = 3.14159265358979d0
      x35 = 1d0/x34
      x36 = 0.009578475d0*x35
      x37 = 0.0676322201645715d0*x12 + 0.00126674656487366d0*x30*x31 - x
     &33*x36 + 0.0188171922990732d0*x4 + 1.0d0
      x38 = x28*x37
      x39 = x27*x38
      x40 = x26*x39
      x41 = rho_c**(-1.66666666666667d0)
      x42 = 0.0223069557172962d0*x41
      x43 = x25*x40*x42
      x44 = x22*x4
      x45 = 5.57236303607474d0*mu*x16 + 1.0d0
      x46 = mu**3
      x47 = 3.39260255780131d0*x14*x46 + 6.76834298980504d0*x44 + x45
      x48 = x47/(3.13317926779378d0*x44 + x45)
      x49 = x1*(x6 - 1.0d0)
      x50 = x49*log(x48)
      x51 = 0.00478594012208448d0*x20
      x52 = rho_c**(-2.66666666666667d0)
      x53 = 0.133841734303777d0*x52
      x54 = mu**8*x51*x53
      x55 = x38 - 1.0d0
      x56 = 0.148394183600699d0
      x57 = x41*x56
      x58 = 0.0102659822546843d0*x20
      x59 = rho_c**(-2.0d0)
      x60 = 1.0525946465532d0*x59
      x61 = mu**6*(0.0167302167879722d0*x55*x57 + x58*x60)
      x62 = 0.179587122125167d0
      x63 = 0.175432441092201d0*x38*x62
      x64 = exp(-0.28165369188898165d0*x4)
      x65 = x64*(0.825481812223657d0*x12 - 4.49737346725955d0*x4)
      x66 = 12.0d0*x34
      x67 = x65*x66
      x68 = 0.0524148278841779d0*x12 + 0.494402081358784d0*x4 + 1.0d0
      x69 = 1d0/x68
      x70 = -0.0259335011650457d0*x4 + 1.0d0
      x71 = x69*x70
      x72 = 2.28942848510666d0
      x73 = x29*x72
      x74 = 1.63540853354893d0*x73
      x75 = x71*x74
      x76 = x26*(x10*x75 + x32*x67)
      x77 = 0.00111534778586481d0*x41
      x78 = x27*x46*(x33*x63 + x76*x77)
      x79 = mu**4
      x80 = 0.13157433081915d0*x55
      x81 = x33*x35
      x82 = 0.022020833726518d0*x20
      x83 = 3.10428814221102d0*x31
      x84 = exp(-0.49698248213959023d0*x4)
      x85 = x84*(0.558025705063192d0*x12 - 0.352521395009435d0*x4)
      x86 = x8*x85
      x87 = x56*(1.63540853354893d0*x73 - x75 - 10.9027235569928d0*x86)
      x88 = x79*(0.00139418473233101d0*rho_c**(-1.0d0)*x87 - x80*x81 - x
     &82*x83)
      x89 = x24*(x43 - x50 + x54 + x61 + x78 - x88)
      x90 = x41*x9
      x91 = rho_c**(-1.5d0)*x13
      x92 = x3*x31
      x93 = rho_c**(-1.16666666666667d0)*x15
      x94 = x5*(0.27127533634502d0*x90 + 0.709361408239834d0*x91 + 1.086
     &51697314076d0*x92 + 1.20668365571947d0*x93)/(x17**2*(x18 + 2.0d0))
     &
      x95 = rho_c**(-2.33333333333333d0)
      x96 = x30*x95
      x97 = x28*(x36*x59 - 0.045088146776381d0*x90 - 0.00627239743302441
     &d0*x92 - 0.00168899541983155d0*x96)
      x98 = rho_c**(-3.0d0)
      x99 = mu**7
      x100 = 0.928727172679123d0*x93
      x101 = mu*x92
      x102 = x35*x59
      x103 = x28*(0.0047892375d0*x102 + 0.113935126864452d0*x37*x92 - 0.
     &0225440733881905d0*x90 - 0.0031361987165122d0*x92 - 0.000844497709
     &915773d0*x96)
      x104 = 0.0141372897033723d0*x11*x34*x69*x72
      x105 = rho_c**(-0.333333333333333d0)
      x106 = x105*x73
      x107 = 1.09027235569928d0*x106*x71
      x108 = x10*x70*x74*(0.034943218589452d0*x90 + 0.164800693786261d0*
     &x92)/x68**2
      E = -rho_c*(x21 - x89)
      d1E(1) = -rho_c*(mu*x24*(2.83060463816878d-5*rho_c**(-4.0d0)*x19*x
     &99 + 0.356911291476739d0*rho_c**(-3.66666666666667d0)*x51*x99 + x2
     &2* x27*(-x26*x77*(2.41662179562687d0*rho_c**(-0.333333333333333d0)
     &* x65 + x104 + x107 + x108 + x32*x64*x66*(-0.550321208149104d0*x90
     & + 1.49912448908652d0*x92) + x67) - 0.175432441092201d0*x33*x62* x
     &97 - 0.00490180588786583d0*x38*x95 + 0.00185891297644135d0*x52* x7
     &6 + x59*x63) + x25*(0.00047750955226877d0*rho_c**( -3.333333333333
     &33d0)*x19 - 0.0334604335759443d0*x103*x57 + 0.0278836946466203d0*x
     &52*x55*x56 + 2.10518929310641d0*x58*x98 - 0.101321183642338d0*x60*
     &x94) - x26*x27*x42*x79*x97 - 0.000290571663240495d0*x39*x79*x98 + 
     &0.037178259528827d0*x40*x52* x79 + x46*(x102*x80 - 0.2631486616383
     &01d0*x103*x81 + 0.00302075971817057d0*x19*x52 + 0.0013941847323310
     &1d0*x57*( -10.9027235569928d0*x10*x8*x84*(-0.372017136708795d0*x90
     & + 0.117507131669812d0*x92) - x104 - 7.26848237132856d0*x105*x86 +
     & 1.09027235569928d0*x106 - x107 - x108 - 1.80615420514536d0*x11*x2
     & *x85) - 0.00232364122055169d0*x59*x87 + 4.13905085628136d0*x82* x
     &95 - 0.217336917462899d0*x83*x94) - 0.0472353356922751d0*x53*x94 *
     &x99 + x49*(-x100 - 2.25611432993501d0*x101 - 1.69630127890066d0* x
     &22*x91 + x48*(x100 + 1.04439308926459d0*x101))/x47) - 0.0089546919
     &0170509d0*x19*x31 + 1.35631049481572d0*x22*x90*(-x43 + x50 - x54 -
     & x61 - x78 + x88)/x23**5 + 2.0d0*x94) - x21 + x89
      end subroutine


C*****************************************************************************
      pure subroutine D2ESRC_PW92_ERF(rho_c, mu, E, d1E, d2E)
C*****************************************************************************
C   Implemented by E.R. Kjellgren.
C
C   Subroutine generated using Sympy 1.3
C   Generated: March 25, 2019
C*****************************************************************************
      implicit none
      real*8, intent(in) :: rho_c, mu
      real*8, intent(out) :: E, d1E(9), d2E(45)
      real*8 :: x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x
     &13, x14, x15, x16, x17, x18, x19, x20, x21, x22, x23, x24, x25, x2
     &6, x27, x28, x29, x30, x31, x32, x33, x34, x35, x36, x37, x38, x39
     &, x40, x41, x42, x43, x44, x45, x46, x47, x48, x49, x50, x51, x52,
     & x53, x54, x55, x56, x57, x58, x59, x60, x61, x62, x63, x64, x65, 
     &x66, x67, x68, x69, x70, x71, x72, x73, x74, x75, x76, x77, x78, x
     &79, x80, x81, x82, x83, x84, x85, x86, x87, x88, x89, x90, x91, x9
     &2, x93, x94, x95, x96, x97, x98, x99, x100, x101, x102, x103, x104
     &, x105, x106, x107, x108, x109, x110, x111, x112, x113, x114, x115
     &, x116, x117, x118, x119, x120, x121, x122, x123, x124, x125, x126
     &, x127, x128, x129, x130, x131, x132, x133, x134, x135, x136, x137
     &, x138, x139, x140, x141, x142, x143, x144, x145, x146, x147, x148
     &, x149, x150, x151, x152, x153, x154, x155, x156, x157, x158, x159
     &, x160, x161, x162, x163, x164, x165, x166, x167, x168, x169, x170
     &, x171, x172, x173, x174, x175, x176, x177, x178, x179, x180, x181
     &, x182, x183, x184, x185, x186, x187, x188, x189, x190, x191, x192
     &, x193, x194, x195, x196, x197, x198, x199, x200, x201, x202, x203
     &, x204, x205, x206, x207, x208, x209, x210, x211, x212, x213, x214
     &, x215, x216, x217, x218, x219, x220, x221, x222, x223, x224, x225
     &, x226, x227, x228, x229, x230, x231, x232, x233, x234, x235, x236
     &, x237, x238, x239, x240, x241, x242, x243, x244
      E = 0.0d0
      d1E(:) = 0.0d0
      d2E(:) = 0.0d0
      x0 = 9.86960440108936d0
      x1 = 2.0d0/x0
      x2 = 1.46459188756152d0
      x3 = 1d0/x2
      x4 = rho_c**(-0.333333333333333d0)*x3
      x5 = 0.194159335344114d0*x4 + 1.0d0
      x6 = 0.693147180559945d0
      x7 = -x6 + 1.0d0
      x8 = 2.14502939711103d0
      x9 = 1d0/x8
      x10 = rho_c**0.666666666666667d0
      x11 = 1d0/x10
      x12 = x11*x9
      x13 = rho_c**(-0.5d0)
      x14 = 0.564189583547756d0
      x15 = 1.41872281647967d0*x14
      x16 = 0.826307487110758d0
      x17 = rho_c**(-0.166666666666667d0)*x16
      x18 = 0.406913004517529d0*x12 + x13*x15 + 7.24010193431683d0*x17 +
     & 3.25955091942229d0*x4
      x19 = 1d0/x7
      x20 = x0*x19
      x21 = x20/x18
      x22 = x7*log(0.5d0*x21 + 1.0d0)
      x23 = x22*x5
      x24 = x1*x23
      x25 = mu**2
      x26 = 0.508616435555896d0*x12*x25 + 1.0d0
      x27 = x26**(-4)
      x28 = mu**5
      x29 = 1.41421356237310d0
      x30 = 0.0837224526465878d0
      x31 = exp(-0.6836107611867116d0*x4)
      x32 = 4.60115111447049d0
      x33 = 1d0/x32
      x34 = rho_c**(-1.33333333333333d0)
      x35 = rho_c**1.0d0
      x36 = 1d0/x35
      x37 = 3.14159265358979d0
      x38 = 1d0/x37
      x39 = 0.009578475d0*x38
      x40 = 0.0676322201645715d0*x12 + 0.00126674656487366d0*x33*x34 - x
     &36*x39 + 0.0188171922990732d0*x4 + 1.0d0
      x41 = x31*x40
      x42 = x30*x41
      x43 = x29*x42
      x44 = rho_c**(-1.66666666666667d0)
      x45 = 0.0223069557172962d0*x44
      x46 = x28*x43*x45
      x47 = x25*x4
      x48 = 5.57236303607474d0*mu*x17 + 1.0d0
      x49 = 3.13317926779378d0*x47 + x48
      x50 = 1d0/x49
      x51 = mu**3
      x52 = 3.39260255780131d0*x13*x14*x51 + 6.76834298980504d0*x47 + x4
     &8
      x53 = x50*x52
      x54 = x1*(x6 - 1.0d0)
      x55 = x54*log(x53)
      x56 = 0.00478594012208448d0*x23
      x57 = rho_c**(-2.66666666666667d0)
      x58 = 0.133841734303777d0*x57
      x59 = mu**8*x56*x58
      x60 = x41 - 1.0d0
      x61 = 0.148394183600699d0
      x62 = x44*x61
      x63 = 0.0102659822546843d0*x23
      x64 = rho_c**(-2.0d0)
      x65 = 1.0525946465532d0*x64
      x66 = mu**6*(0.0167302167879722d0*x60*x62 + x63*x65)
      x67 = 0.179587122125167d0
      x68 = x41*x67
      x69 = 0.175432441092201d0*x68
      x70 = 0.825481812223657d0*x12 - 4.49737346725955d0*x4
      x71 = exp(-0.28165369188898165d0*x4)
      x72 = x37*x71
      x73 = 12.0d0*x72
      x74 = x70*x73
      x75 = 0.0524148278841779d0*x12 + 0.494402081358784d0*x4 + 1.0d0
      x76 = 1d0/x75
      x77 = -0.0259335011650457d0*x4 + 1.0d0
      x78 = x76*x77
      x79 = 2.28942848510666d0
      x80 = x32*x79
      x81 = 1.63540853354893d0*x80
      x82 = x78*x81
      x83 = x30*(x10*x82 + x35*x74)
      x84 = 0.00111534778586481d0*x44
      x85 = x29*x51*(x36*x69 + x83*x84)
      x86 = mu**4
      x87 = 0.13157433081915d0*x60
      x88 = x36*x38
      x89 = 0.022020833726518d0*x23
      x90 = 3.10428814221102d0*x34
      x91 = 1.63540853354893d0*x80
      x92 = exp(-0.49698248213959023d0*x4)
      x93 = x92*(0.558025705063192d0*x12 - 0.352521395009435d0*x4)
      x94 = x8*x93
      x95 = 10.9027235569928d0*x94
      x96 = x61*(-x82 + x91 - x95)
      x97 = x86*(0.00139418473233101d0*rho_c**(-1.0d0)*x96 - x87*x88 - x
     &89*x90)
      x98 = x46 - x55 + x59 + x66 + x85 - x97
      x99 = x27*x98
      x100 = 0.0691804894611514d0
      x101 = x100*x22
      x102 = 0.129439556896076d0*x34
      x103 = 2.0d0*x5
      x104 = x21 + 2.0d0
      x105 = 1d0/x104
      x106 = x18**(-2)
      x107 = x44*x9
      x108 = rho_c**(-1.5d0)
      x109 = x108*x14
      x110 = x3*x34
      x111 = rho_c**(-1.16666666666667d0)*x16
      x112 = 0.27127533634502d0*x107 + 0.709361408239834d0*x109 + 1.0865
     &1697314076d0* x110 + 1.20668365571947d0*x111
      x113 = x105*x106*x112
      x114 = x26**(-5)
      x115 = x107*x114
      x116 = -x46 + x55 - x59 - x66 - x85 + x97
      x117 = x116*x25
      x118 = rho_c**(-0.333333333333333d0)
      x119 = x118*x80
      x120 = 1.09027235569928d0*x119
      x121 = x11*x2
      x122 = 1.80615420514536d0*x121*x93
      x123 = 7.26848237132856d0*x118*x94
      x124 = x92*(-0.372017136708795d0*x107 + 0.117507131669812d0*x110)
      x125 = x124*x8
      x126 = 10.9027235569928d0*x10
      x127 = x125*x126
      x128 = x11*x37*x79
      x129 = 0.0141372897033723d0*x128*x76
      x130 = 1.09027235569928d0*x119*x78
      x131 = x75**(-2)
      x132 = 0.034943218589452d0*x107 + 0.164800693786261d0*x110
      x133 = x131*x132
      x134 = x133*x77
      x135 = x10*x81
      x136 = x134*x135
      x137 = x120 - x122 - x123 - x127 - x129 - x130 - x136
      x138 = 0.00139418473233101d0*x62
      x139 = 0.00232364122055169d0*x64
      x140 = x38*x64
      x141 = 0.0150354743280612d0*x22
      x142 = 0.200908840802769d0*x57
      x143 = rho_c**(-2.33333333333333d0)
      x144 = x143*x33
      x145 = -0.0225440733881905d0*x107 - 0.0031361987165122d0*x110 + 0.
     &0047892375d0* x140 - 0.000844497709915773d0*x144
      x146 = x31*(0.113935126864452d0*x110*x40 + x145)
      x147 = 0.263148661638301d0*x88
      x148 = 0.217336917462899d0
      x149 = x113*x5
      x150 = x148*x149
      x151 = x140*x87 + x141*x142 + 4.13905085628136d0*x143*x89 - x146*x
     &147 - x150* x90
      x152 = 0.00168899541983155d0*x144
      x153 = x39*x64
      x154 = 0.045088146776381d0*x107
      x155 = 0.00627239743302441d0*x110
      x156 = x31*(-x152 + x153 - x154 - x155)
      x157 = x156*x67
      x158 = 0.175432441092201d0*x36
      x159 = x157*x158
      x160 = 0.122619224952946d0
      x161 = x160*x41
      x162 = 0.0399758348639607d0*x143*x161
      x163 = x64*x69
      x164 = 0.00185891297644135d0*x57*x83
      x165 = -0.550321208149104d0*x107 + 1.49912448908652d0*x110
      x166 = x35*x73
      x167 = x70*x71
      x168 = 2.14502939711103d0
      x169 = rho_c**(-0.333333333333333d0)*x168
      x170 = x129 + x130 + x136
      x171 = x30*(x165*x166 + 1.12661476755593d0*x167*x169 + x170 + x74)
     &
      x172 = x171*x84
      x173 = x25*x29
      x174 = 0.928727172679123d0*x111
      x175 = mu*x110
      x176 = 2.25611432993501d0*x175
      x177 = 1.69630127890066d0*x109*x25
      x178 = x174 + 1.04439308926459d0*x175
      x179 = -x174 - x176 - x177 + x178*x53
      x180 = 1d0/x52
      x181 = x180*x54
      x182 = x179*x181
      x183 = x156*x30
      x184 = x29*x86
      x185 = x184*x45
      x186 = x57*x61
      x187 = 0.00700944907716014d0*x22
      x188 = rho_c**(-3.33333333333333d0)
      x189 = 0.0681236923205142d0*x188
      x190 = rho_c**(-3.0d0)
      x191 = 0.0334604335759443d0*x62
      x192 = 0.101321183642338d0
      x193 = x192*x65
      x194 = mu**7
      x195 = 0.00326776364305339d0*x194*x22
      x196 = rho_c**(-4.0d0)
      x197 = 0.00866220739124163d0*x196
      x198 = rho_c**(-3.66666666666667d0)
      x199 = x194*x56
      x200 = x57*x86
      x201 = 0.0571643564037363d0
      x202 = x184*x41
      x203 = x201*x202
      x204 = x194*x58
      x205 = 0.0472353356922751d0
      x206 = x149*x205
      x207 = -0.00508309165921971d0*x190*x203 + x195*x197 + 0.3569112914
     &76739d0*x198* x199 + 0.037178259528827d0*x200*x43 - x204*x206 + x2
     &8*(-x146*x191 - x149*x193 + 0.0278836946466203d0*x186*x60 + x187*x
     &189 + 2.10518929310641d0*x190*x63)
      x208 = x173*(-x159 - x162 + x163 + x164 - x172) + x182 - x183*x185
     & + x207 + x51 *(x137*x138 - x139*x96 + x151)
      x209 = mu*x27
      x210 = 2.71262098963145d0*x115
      x211 = x185*x30
      x212 = x103*x105
      x213 = x57*x9
      x214 = rho_c**(-2.5d0)*x14
      x215 = x143*x3
      x216 = rho_c**(-2.16666666666667d0)*x16
      x217 = x106*(0.452125560575033d0*x213 + 1.06404211235975d0*x214 + 
     &1.44868929752102d0*x215 + 1.40779759833938d0*x216)
      x218 = x112**2/(x104**2*x18**4)
      x219 = x112*(0.542550672690039d0*x107 + x108*x15 + 2.1730339462815
     &3d0*x110 + 2.41336731143894d0*x111)/x18**3
      x220 = x188*x33
      x221 = x190*x38
      x222 = x31*(0.0751469112939684d0*x213 + 0.00836319657736588d0*x215
     & + 0.00394098931294027d0*x220 - 0.01915695d0*x221)
      x223 = rho_c**(-4.33333333333333d0)
      x224 = x113*x194*x197
      x225 = x105*x217
      x226 = x204*x205*x5
      x227 = x19*x218
      x228 = x227*x5
      x229 = x105*x219
      x230 = x174 + x176 + x177
      x231 = mu*x178
      x232 = x231*x50
      x233 = 1.08351503479231d0*x216
      x234 = mu*x215
      x235 = rho_c**(-3.66666666666667d0)
      x236 = rho_c**(-1.33333333333333d0)
      x237 = x236*x80
      x238 = x135*x77
      x239 = -2.18054471139857d0*x119*x134 - 0.0282745794067445d0*x128*x
     &133 + x131* x238*(0.0582386976490866d0*x213 + 0.219734258381682d0*
     &x215) - x132*x238*(0.0698864371789039d0*x107 + 0.329601387572523d0
     &*x110)/ x75**3 + 0.363424118566428d0*x237*x78
      x240 = x113*x189
      x241 = x193*x5
      x242 = x31*(0.455740507457808d0*x110*x145 - x190*x39 + 0.025962426
     &2672375d0* x213*x40 + 0.0375734556469842d0*x213 - 0.15191350248593
     &6d0*x215* x40 + 0.00418159828868294d0*x215 + 0.00197049465647014d0
     &*x220)
      x243 = x5*x90
      x244 = x148*x243
      E = -rho_c*(x24 - x99)
      d1E(1) = -rho_c*(-x101*x102 + x103*x113 + 1.35631049481572d0*x115*
     &x117 + x208* x209) - x24 + x99
      d2E(1) = -rho_c*(0.172586075861435d0*x101*x143 - 0.682784063255296
     &d0*x102*x113 - x103*x20*x218 - 0.129439556896076d0*x110*x113 - 2.2
     &6051749135954d0*x114*x117*x213 + 2.29947269793409d0*x116*x220* x86
     &/x26**6 + x208*x210*x51 + x209*(mu*x179*x230*x54/x52**2 - 0.057748
     &0492749442d0*rho_c**(-5.0d0)*x195 - 1.30867473541471d0* rho_c**(-4
     &.66666666666667d0)*x199 - 0.0101661833184394d0*x156* x184*x190*x20
     &1 + x173*(-0.0799516697279215d0*x143*x156*x160 + 0.350864882184401
     &d0*x157*x64 - x158*x222*x67 + 0.133252782879869d0*x161*x188 + 0.00
     &37178259528827d0*x171*x57 - 0.350864882184401d0*x190*x68 - 0.00495
     &710127051027d0*x198*x83 - 0.00910930363347549d0*x235*x42 + x30*x84
     &*(-2.25322953511185d0* x165*x169*x71 - 24.0d0*x165*x72 - x166*(0.9
     &17202013581841d0*x213 - 1.99883265211535d0*x215) - 0.7510765117039
     &51d0*x167*x168*x34 - 0.154912426780983d0*x167*x44 + x239)) + x181*
     &(2.54445191835098d0* x214*x25 - 2.0d0*x230*x232 + x231*x52*(1.8574
     &5434535825d0*x111 + 2.08878617852919d0*x175)/x49**2 + x233 + 3.008
     &15243991335d0*x234 - x53*(x233 + 1.39252411901946d0*x234)) - x182*
     &x232 + 0.074356519057654d0*x183*x200*x29 + 0.713822582953479d0*x19
     &4*x198 *x206 + 0.023721094409692d0*x196*x203 - 0.0991420254102054d
     &0*x198 *x43*x86 - 4.52089344419912d-5*x202*x223 + 0.46619407703541
     &2d0* x204*x228 - x211*x222 + 0.064503068866399d0*x224 + x225*x226 
     &- x226*x229 + x28*(x100*x240 + 0.111534778586481d0*x146*x186 + 4.2
     &1037858621281d0*x149*x190*x192 - 0.363326359042743d0*x187*x223 - x
     &191*x242 - 6.31556787931922d0*x196*x63 - 0.074356519057654d0* x198
     &*x60*x61 + x225*x241 + x228*x65 - x229*x241 + 0.0691804894611514d0
     &*x240) + x51*(0.148394183600699d0*x113*x142 + 0.200908840802769d0*
     &x113*x186 - 0.00464728244110338d0*x137*x186 + x138*(-14.5369647426
     &571d0*x118*x125 - 3.61230841029072d0*x121* x124 - x126*x8*x92*(0.6
     &20028561181324d0*x213 - 0.156676175559749d0*x215) + 2.422827457109
     &52d0*x236*x94 - 0.363424118566428d0*x237 + x239 - 0.299209d0*x64*x
     &93) + 0.526297323276602d0*x140*x146 - 0.535756908807384d0*x141*x19
     &8 - 0.267878454403692d0*x141*x235 + 8.27810171256273d0*x143*x150 -
     & x147*x242 - 9.65778533132318d0*x188*x89 + 0.00619637658813783d0* 
     &x190*x96 - 0.263148661638301d0*x221*x60 + x225*x244 + 2.1450293971
     &1103d0*x227*x243 - x229*x244)) - x212*x217 + x212* x219) + 0.25887
     &9113792152d0*x101*x34 - 4.0d0*x149 - 2.0d0*x209*( -x1*x179*x180*x7
     & - x173*(x159 + x162 - x163 - x164 + x172) + x207 + x211*x31*(x152
     & - x153 + x154 + x155) + x51*(-x138*(-x120 + x122 + x123 + x127 + 
     &x170) + x139*x61*(x82 - x91 + x95) + x151)) + x210*x25*x98
      end subroutine


C*****************************************************************************
      pure subroutine ESRC_SPIN_PW92_ERF(rho_c, rho_s, mu, E)
C*****************************************************************************
C   Implemented by E.R. Kjellgren.
C
C   Subroutine generated using Sympy 1.3
C   Generated: March 25, 2019
C*****************************************************************************
      implicit none
      real*8, intent(in) :: rho_c, rho_s, mu
      real*8, intent(out) :: E
      real*8 :: x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x
     &13, x14, x15, x16, x17, x18, x19, x20, x21, x22, x23, x24, x25, x2
     &6, x27, x28, x29, x30, x31, x32, x33, x34, x35, x36, x37, x38, x39
     &, x40, x41, x42, x43, x44, x45, x46, x47, x48, x49, x50
      E = 0.0d0
      x0 = 0.682784063255296d0*rho_c**(-0.333333333333333d0)
      x1 = 2.14502939711103d0
      x2 = rho_c**0.666666666666667d0
      x3 = 1/(x1*x2)
      x4 = 0.564189583547756d0*rho_c**(-0.5d0)
      x5 = 0.826307487110758d0*rho_c**(-0.166666666666667d0)
      x6 = 9.86960440108936d0
      x7 = 0.306852819440055d0
      x8 = x6/x7
      x9 = 2.0d0*(0.194159335344114d0*x0 + 1.0d0)*log(0.5d0*x8/( 3.25955
     &09194222925d0*x0 + 0.40691300451752932d0*x3 + 1.4187228164796675d0
     &*x4 + 7.2401019343168302d0*x5) + 1.0d0)
      x10 = x7/x6
      x11 = x10*x9
      x12 = rho_s**4/rho_c**4
      x13 = rho_s/rho_c
      x14 = x13 + 1.0d0
      x15 = -x13 + 1.0d0
      x16 = 1.92366105093154d0*x14**1.33333333333333d0 + 1.9236610509315
     &4d0*x15** 1.33333333333333d0 - 3.84732210186307d0
      x17 = 0.0197517897025652d0*x16*(0.101077332976288d0*x0 + 1.0d0)*(-
     &x12 + 1.0d0) *log(1.0d0 + 29.608574643216677d0/(3.2918048099450625
     &d0*x0 + 0.41002507094961249d0*x3 + 0.76232752193529008d0*x4 + 9.87
     &21297225692712d0*x5))
      x18 = x10*x12*x16*(-x9 + (0.186690969707574d0*x0 + 1.0d0)*log(x8/(
     & 5.630984149097876d0*x0 + 0.51606646454786342d0*x3 + 2.91521471421
     &91778d0*x4 + 13.457913714394447d0*x5) + 1.0d0))
      x19 = mu**2
      x20 = 1.41421356237310d0
      x21 = rho_c**(-1.66666666666667d0)
      x22 = 0.0837224526465878d0*x21
      x23 = 1.0d0 - rho_s**2/rho_c**2
      x24 = 4.60115111447049d0
      x25 = rho_c**(-1.33333333333333d0)/x24
      x26 = 3.14159265358979d0
      x27 = rho_c**1.0d0
      x28 = 1d0/x27
      x29 = x28/x26
      x30 = (0.0188171922990732d0*x0 + 0.00126674656487366d0*x25 - 0.009
     &578475d0*x29 + 0.0676322201645715d0*x3 + 1.0d0)*exp(-0.68361076118
     &67116d0*x0)
      x31 = x23*x30
      x32 = x14**0.666666666666667d0 + x15**0.666666666666667d0
      x33 = x32**3
      x34 = x0*x19/x32**2
      x35 = 11.1447260721495d0*mu*x5/x32 + 1.0d0
      x36 = mu**3
      x37 = 1d0/x14
      x38 = x37**0.666666666666667d0
      x39 = 0.0524148278841779d0*x3
      x40 = x0*x37**0.333333333333333d0
      x41 = 2.28942848510666d0*x24
      x42 = 0.817704266774463d0*x41
      x43 = x14**2*x42*(-0.0259335011650457d0*x40 + 1.0d0)/(x38*(x38*x39
     & + 0.494402081358784d0*x40 + 1.0d0))
      x44 = 1d0/x15
      x45 = x44**0.666666666666667d0
      x46 = x0*x44**0.333333333333333d0
      x47 = x15**2*x42*(-0.0259335011650457d0*x46 + 1.0d0)/(x45*(x39*x45
     & + 0.494402081358784d0*x46 + 1.0d0))
      x48 = -x11 + x17 - x18
      x49 = 0.148394183600699d0
      x50 = x23*(x30 - 1.0d0)
      E = -rho_c*(x11 - x17 + x18 + 0.25d0*(0.0126441184989504d0*mu**8*r
     &ho_c**( -2.66666666666667d0)*x48 + 4.0d0*mu**6*(0.0533250677421794
     &d0* rho_c**(-2.0d0)*x48 - 0.0167302167879722d0*x21*x49*x50) - 0.08
     &92278228691848d0*mu**5*x20*x22*x31 + 4.0d0*mu**4*( -0.001394184732
     &33101d0*rho_c**(-1.0d0)*x49*(10.9027235569928d0*x1 *x23*(-0.352521
     &395009435d0*x0 + 0.558025705063192d0*x3)*exp( -0.49698248213959023
     &d0*x0) - 0.817704266774463d0*x41*(x14** 2.66666666666667d0 + x15**
     &2.66666666666667d0) + x43 + x47) + 1.55214407110551d0*x25*x48 - 0.
     &13157433081915d0*x29*x50) - x10* x33*log((27.073371959220147d0*x34
     & + x35 + 27.140820462410485d0* x36*x4/x33)/(12.532717071175123d0*x
     &34 + x35)) - 4.0d0*x20*x36*( 0.00111534778586481d0*x22*(x2*x43 + x
     &2*x47 + 12.0d0*x23*x26*x27*( -4.49737346725955d0*x0 + 0.8254818122
     &23657d0*x3)*exp( -0.28165369188898165d0*x0)) + 0.0315054072231411d
     &0*x28*x31))/( 0.508616435555896d0*x19*x3 + 1.0d0)**4)
      end subroutine


C*****************************************************************************
      pure subroutine D1ESRC_SPIN_PW92_ERF(rho_c, rho_s, mu, E, d1E)
C*****************************************************************************
C   Implemented by E.R. Kjellgren.
C
C   Subroutine generated using Sympy 1.3
C   Generated: March 25, 2019
C*****************************************************************************
      implicit none
      real*8, intent(in) :: rho_c, rho_s, mu
      real*8, intent(out) :: E, d1E(9)
      real*8 :: x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x
     &13, x14, x15, x16, x17, x18, x19, x20, x21, x22, x23, x24, x25, x2
     &6, x27, x28, x29, x30, x31, x32, x33, x34, x35, x36, x37, x38, x39
     &, x40, x41, x42, x43, x44, x45, x46, x47, x48, x49, x50, x51, x52,
     & x53, x54, x55, x56, x57, x58, x59, x60, x61, x62, x63, x64, x65, 
     &x66, x67, x68, x69, x70, x71, x72, x73, x74, x75, x76, x77, x78, x
     &79, x80, x81, x82, x83, x84, x85, x86, x87, x88, x89, x90, x91, x9
     &2, x93, x94, x95, x96, x97, x98, x99, x100, x101, x102, x103, x104
     &, x105, x106, x107, x108, x109, x110, x111, x112, x113, x114, x115
     &, x116, x117, x118, x119, x120, x121, x122, x123, x124, x125, x126
     &, x127, x128, x129, x130, x131, x132, x133, x134, x135, x136, x137
     &, x138, x139, x140, x141, x142, x143, x144, x145, x146, x147, x148
     &, x149, x150, x151, x152, x153, x154, x155, x156, x157, x158, x159
     &, x160, x161, x162, x163, x164, x165, x166, x167, x168, x169, x170
     &, x171, x172, x173, x174, x175, x176, x177, x178, x179, x180, x181
     &, x182, x183, x184, x185, x186, x187, x188, x189, x190, x191, x192
     &, x193, x194, x195, x196, x197, x198, x199, x200, x201, x202, x203
     &, x204, x205, x206, x207, x208, x209, x210, x211, x212, x213, x214
     &, x215, x216, x217, x218, x219, x220, x221, x222, x223, x224, x225
     &, x226, x227, x228, x229, x230, x231, x232, x233, x234, x235, x236
     &, x237, x238, x239, x240, x241, x242, x243, x244, x245, x246, x247
     &, x248, x249, x250, x251, x252, x253, x254, x255, x256, x257, x258
     &, x259, x260, x261, x262, x263, x264, x265, x266, x267, x268, x269
     &, x270, x271, x272
      E = 0.0d0
      d1E(:) = 0.0d0
      x0 = 2.14502939711103d0
      x1 = 1d0/x0
      x2 = rho_c**0.666666666666667d0
      x3 = 1d0/x2
      x4 = x1*x3
      x5 = 0.564189583547756d0
      x6 = rho_c**(-0.5d0)*x5
      x7 = 1.46459188756152d0
      x8 = 1d0/x7
      x9 = rho_c**(-0.333333333333333d0)*x8
      x10 = 0.826307487110758d0
      x11 = rho_c**(-0.166666666666667d0)*x10
      x12 = 7.24010193431683d0*x11 + 0.406913004517529d0*x4 + 1.41872281
     &647967d0*x6 + 3.25955091942229d0*x9
      x13 = 9.86960440108936d0
      x14 = 0.306852819440055d0
      x15 = x13/x14
      x16 = x15/x12
      x17 = log(0.5d0*x16 + 1.0d0)
      x18 = 0.194159335344114d0*x9 + 1.0d0
      x19 = 2.0d0*x18
      x20 = x17*x19
      x21 = x14/x13
      x22 = x20*x21
      x23 = 1d0/rho_c
      x24 = rho_s*x23
      x25 = x24 + 1.0d0
      x26 = -x24 + 1.0d0
      x27 = 1.92366105093154d0*x25**1.33333333333333d0 + 1.9236610509315
     &4d0*x26** 1.33333333333333d0 - 3.84732210186307d0
      x28 = rho_c**(-4)
      x29 = rho_s**4
      x30 = x28*x29
      x31 = -x30 + 1.0d0
      x32 = 0.101077332976288d0*x9 + 1.0d0
      x33 = 9.87212972256927d0*x11 + 0.410025070949612d0*x4 + 0.76232752
     &193529d0*x6 + 3.29180480994506d0*x9
      x34 = 1.0d0 + 29.6085746432167d0/x33
      x35 = log(x34)
      x36 = x31*x32*x35
      x37 = 0.0197517897025652d0*x27*x36
      x38 = 0.186690969707574d0*x9 + 1.0d0
      x39 = 13.4579137143944d0*x11 + 0.516066464547863d0*x4 + 2.91521471
     &421918d0*x6 + 5.63098414909788d0*x9
      x40 = x15/x39 + 1.0d0
      x41 = log(x40)
      x42 = x38*x41
      x43 = -x20 + x42
      x44 = x27*x30
      x45 = x21*x44
      x46 = x43*x45
      x47 = rho_c**(-2)
      x48 = rho_s**2
      x49 = -x47*x48 + 1.0d0
      x50 = exp(-0.6836107611867116d0*x9)
      x51 = 4.60115111447049d0
      x52 = 1d0/x51
      x53 = rho_c**(-1.33333333333333d0)
      x54 = x52*x53
      x55 = rho_c**1.0d0
      x56 = 1d0/x55
      x57 = 3.14159265358979d0
      x58 = 1d0/x57
      x59 = 0.009578475d0*x58
      x60 = 0.0676322201645715d0*x4 + 0.00126674656487366d0*x54 - x56*x5
     &9 + 0.0188171922990732d0*x9 + 1.0d0
      x61 = x50*x60
      x62 = x49*x61
      x63 = 0.0837224526465878d0
      x64 = 1.41421356237310d0
      x65 = mu**5
      x66 = x63*x64*x65
      x67 = x62*x66
      x68 = rho_c**(-1.66666666666667d0)
      x69 = 0.0892278228691848d0*x68
      x70 = x25**0.666666666666667d0 + x26**0.666666666666667d0
      x71 = x70**3
      x72 = mu**2
      x73 = x70**2
      x74 = 1d0/x73
      x75 = x72*x74
      x76 = x75*x9
      x77 = 1d0/x70
      x78 = mu*x77
      x79 = 11.1447260721495d0*x11*x78 + 1.0d0
      x80 = mu**3
      x81 = 1d0/x71
      x82 = 27.1408204624105d0*x6*x80*x81 + 27.0733719592201d0*x76 + x79
     &
      x83 = x82/(12.5327170711751d0*x76 + x79)
      x84 = x21*log(x83)
      x85 = 0.179587122125167d0
      x86 = 0.175432441092201d0*x56*x85
      x87 = 0.825481812223657d0*x4 - 4.49737346725955d0*x9
      x88 = exp(-0.28165369188898165d0*x9)
      x89 = x57*x87*x88
      x90 = 12.0d0*x49
      x91 = x89*x90
      x92 = 2.28942848510666d0
      x93 = x51*x92
      x94 = 0.817704266774463d0*x93
      x95 = x25**2
      x96 = 1d0/x25
      x97 = x96**0.666666666666667d0
      x98 = 0.0524148278841779d0*x4
      x99 = x96**0.333333333333333d0
      x100 = x9*x99
      x101 = 0.494402081358784d0*x100 + x97*x98 + 1.0d0
      x102 = 1d0/x101
      x103 = x102*x95
      x104 = 1d0/x97
      x105 = 0.0259335011650457d0*x100
      x106 = x104*(-x105 + 1.0d0)
      x107 = x103*x106
      x108 = x107*x94
      x109 = 1d0/x26
      x110 = x109**0.666666666666667d0
      x111 = 1d0/x110
      x112 = x111*x94
      x113 = x109**0.333333333333333d0
      x114 = x113*x9
      x115 = 0.0259335011650457d0*x114
      x116 = -x115 + 1.0d0
      x117 = x26**2
      x118 = x110*x98 + 0.494402081358784d0*x114 + 1.0d0
      x119 = 1d0/x118
      x120 = x117*x119
      x121 = x116*x120
      x122 = x112*x121
      x123 = x63*(x108*x2 + x122*x2 + x55*x91)
      x124 = 0.00111534778586481d0*x68
      x125 = 4.0d0*x64*x80
      x126 = -x22 + x37
      x127 = x126 - x46
      x128 = 0.0472353356922751d0
      x129 = mu**8
      x130 = x127*x128*x129
      x131 = rho_c**(-2.66666666666667d0)
      x132 = 0.267683468607554d0*x131
      x133 = 0.148394183600699d0
      x134 = x61 - 1.0d0
      x135 = x134*x49
      x136 = x133*x135
      x137 = rho_c**(-2.0d0)
      x138 = 0.101321183642338d0
      x139 = 0.526297323276602d0*x138
      x140 = x137*x139
      x141 = 4.0d0*mu**6
      x142 = 0.13157433081915d0*x135
      x143 = x56*x58
      x144 = x93*(x25**2.66666666666667d0 + x26**2.66666666666667d0)
      x145 = 0.817704266774463d0*x144
      x146 = 0.558025705063192d0*x4 - 0.352521395009435d0*x9
      x147 = exp(-0.49698248213959023d0*x9)
      x148 = x0*x146*x147
      x149 = x148*x49
      x150 = 10.9027235569928d0*x149
      x151 = 0.00139418473233101d0*x133
      x152 = 1.55214407110551d0*x54
      x153 = 4.0d0*mu**4
      x154 = -x125*(x123*x124 + x62*x86) + x130*x132 + x141*(x127*x140 -
     & 0.0167302167879722d0*x136*x68) + x153*(-rho_c**(-1.0d0)*x151*( x1
     &08 + x122 - x145 + x150) + x127*x152 - x142*x143) - x67*x69 - x71*
     &x84
      x155 = 0.508616435555896d0*x4*x72 + 1.0d0
      x156 = x155**(-4)
      x157 = 0.25d0*x156
      x158 = x154*x157
      x159 = 0.0691804894611514d0*x14*x53
      x160 = x159*x17
      x161 = x53*x8
      x162 = x27*x35
      x163 = x161*x162*x31
      x164 = x29/rho_c**5
      x165 = -x25**0.333333333333333d0 + x26**0.333333333333333d0
      x166 = x165*x36
      x167 = rho_s*x47
      x168 = x1*x68
      x169 = rho_c**(-1.5d0)*x5
      x170 = rho_c**(-1.16666666666667d0)*x10
      x171 = (1.08651697314076d0*x161 + 0.27127533634502d0*x168 + 0.7093
     &61408239834d0 *x169 + 1.20668365571947d0*x170)/(x12**2*(x16 + 2.0d
     &0))
      x172 = x27*x31*x32*(1.09726826998169d0*x161 + 0.273350047299742d0*
     &x168 + 0.381163760967645d0*x169 + 1.64535495376155d0*x170)/(x33**2
     &*x34)
      x173 = rho_s**5/rho_c**6
      x174 = x21*x43
      x175 = x164*x27
      x176 = x38*(1.87699471636596d0*x161 + 0.344044309698576d0*x168 + 1
     &.45760735710959d0*x169 + 2.24298561906574d0*x170)/(x39**2*x40)
      x177 = 0.062230323235858d0*x159*x41
      x178 = x171*x19
      x179 = 0.129439556896076d0*x160
      x180 = -x178 + x179
      x181 = x44*(x176 - x177 + x180)
      x182 = rho_c**(-4.66666666666667d0)
      x183 = x48*x61
      x184 = 0.17845564573837d0*x66
      x185 = rho_c**(-2.33333333333333d0)
      x186 = x185*x52
      x187 = x49*x50
      x188 = x187*(-x137*x59 + 0.00627239743302441d0*x161 + 0.0450881467
     &76381d0*x168 + 0.00168899541983155d0*x186)
      x189 = rho_c**(-3.0d0)
      x190 = x25**(-0.333333333333333d0)
      x191 = x26**(-0.333333333333333d0)
      x192 = x190 - x191
      x193 = rho_s*x192
      x194 = 2.0d0*x73*x84
      x195 = rho_c**(-3.66666666666667d0)
      x196 = 1.85745434535825d0*x170
      x197 = x161*x78
      x198 = x169*x75
      x199 = 7.42981738143299d0*rho_c**(-2.16666666666667d0)*x10*x193*x7
     &7
      x200 = rho_c**(-2.33333333333333d0)
      x201 = rho_s*x200
      x202 = mu*x192*x201*x74*x8
      x203 = mu*x21/x82
      x204 = 4.0d0*x175
      x205 = 2.56488140124205d0*x165
      x206 = x173*x205
      x207 = 0.0790071588102608d0*x162*x32
      x208 = 0.0506609980493537d0*x166
      x209 = -0.000665486074881265d0*x163 + x164*x207 + x167*x208 + 0.58
     &482233974552d0*x172 + x180
      x210 = x128*x129
      x211 = rho_c**(-4.0d0)
      x212 = 0.350864882184401d0*x85
      x213 = 24.0d0*x89
      x214 = 1d0/x99
      x215 = x201*x96
      x216 = 0.00706864485168613d0*x57*x92
      x217 = x2*x216
      x218 = x103*x214*x217*(-x215 + x53)
      x219 = x102*x25
      x220 = rho_c**(-1.33333333333333d0)*rho_s
      x221 = x220*x93
      x222 = 2.18054471139857d0*x221
      x223 = x106*x219*x222
      x224 = rho_c**(-0.333333333333333d0)
      x225 = x224*x93
      x226 = 0.545136177849642d0*x225
      x227 = x107*x226
      x228 = 1d0/x113
      x229 = x109*x201
      x230 = x120*x217*x228*(x229 + x53)
      x231 = x119*x26
      x232 = x111*x231
      x233 = x116*x222*x232
      x234 = x111*x121*x226
      x235 = 0.034943218589452d0*rho_s*x1*x131
      x236 = 0.164800693786261d0*x8
      x237 = 0.034943218589452d0*x168
      x238 = 0.164800693786261d0*x161
      x239 = x237*x97 + x238*x99
      x240 = x2*x94/x101**2
      x241 = x106*x240*x95*(-x215*x236*x99 - x235*x96*x97 + x239)
      x242 = x110*x237 + x113*x238
      x243 = x112*x2/x118**2
      x244 = x116*x117*x243*(x109*x110*x235 + x113*x229*x236 + x242)
      x245 = x124*x63
      x246 = x134*x48
      x247 = 0.0334604335759443d0*x133
      x248 = x137*x58
      x249 = x187*(0.113935126864452d0*x161*x60 - 0.0031361987165122d0*x
     &161 - 0.0225440733881905d0*x168 - 0.000844497709915773d0*x186 + 0.
     &0047892375d0*x248)
      x250 = x20 - x42
      x251 = x126 + x250*x45
      x252 = x21*x250
      x253 = -x204*x252 + x206*x252 + x209 + x44*(-x176 + x177 + x178 - 
     &x179)
      x254 = 0.263148661638301d0*x58
      x255 = 21.8054471139857d0*x148
      x256 = -2.18054471139857d0*x25**1.66666666666667d0 + 2.18054471139
     &857d0*x26** 1.66666666666667d0
      x257 = x147*x49
      x258 = x151*x68
      x259 = rho_s**3
      x260 = x207*x259
      x261 = x174*x205
      x262 = 4.0d0*x174*x259*x27
      x263 = rho_s*x61
      x264 = 7.42981738143299d0*x170
      x265 = rho_c**(-3)
      x266 = -x208 - x260*x265 + x261*x30 - x262*x265
      x267 = rho_s*x134
      x268 = x216*x3
      x269 = 2.18054471139857d0*x225
      x270 = x104*(x105 - 1.0d0)
      x271 = x115 - 1.0d0
      x272 = -x214*x219*x268 + x219*x269*x270 + x228*x231*x268 - x232*x2
     &69*x271 + x239*x240*x25*x270 - x242*x243*x26*x271
      E = -rho_c*(x158 + x22 - x37 + x46)
      d1E(1) = 0.25d0*rho_c*(-1.35631049481572d0*x154*x168*x72/x155**5 -
     & x156*(x125*( 0.00185891297644135d0*x123*x131 + 0.175432441092201d
     &0*x137*x62* x85 - x183*x211*x212 - 0.00490180588786583d0*x185*x62 
     &+ x188*x86 - x245*(2.41662179562687d0*rho_c**(-0.333333333333333d0
     &)*x49*x87* x88 + x137*x213*x48 + x218 - x223 + x227 + x230 + x233 
     &+ x234 + x241 + x244 + x55*x57*x88*x90*(1.49912448908652d0*x161 - 
     &0.550321208149104d0*x168) + x91)) - 0.713822582953479d0*x130*x195 
     &+ 0.148713038115308d0*x131*x67 + x132*x210*(x174*x204 - x174*x206 
     &- x181 + x209) - x141*(-0.0278836946466203d0*x131*x136 + 1.0525946
     &465532d0*x138*x189*x251 - x140*x253 + x182*x246*x247 + x247*x249*x
     &68) - x153*(0.00232364122055169d0*x133*x137*(-x108 - x122 + x145 -
     & x150) - x142*x248 + 0.263148661638301d0*x143*x249 - x152*x253 + 2
     &.06952542814068d0*x186*x251 + x211*x246*x254 - x258* (-10.90272355
     &69928d0*x0*x2*x257*(0.117507131669812d0*x161 - 0.372017136708795d0
     &*x168) + 0.545136177849642d0*x144*x224 - 1.80615420514536d0*x146*x
     &257*x3*x7 - 7.26848237132856d0*x149*x224 - x200*x255*x48 - x218 + 
     &x221*x256 + x223 - x227 - x230 - x233 - x234 - x241 - x244)) - x18
     &2*x183*x184 + x188*x66*x69 - 0.00116228665296198d0*x189*x62*x64*x6
     &5 + x193*x194*x47 - x203*x73 *(54.281640924821d0*rho_c**(-2.5d0)*x
     &193*x5*x72*x81 - x196 - 9.02445731974005d0*x197 - 13.5704102312052
     &d0*x198 + x199 + 36.0978292789602d0*x202 + x83*(x196 + 4.177572357
     &05837d0*x197 - x199 - 16.7102894282335d0*x202))) + 0.5177582275843
     &04d0*x160 + 0.316028635241043d0*x162*x164*x32 - 0.0026619442995250
     &6d0*x163 - 10.2595256049682d0*x165*x173*x174 + 0.202643992197415d0
     &*x166*x167 - 8.0d0*x171*x18 + 2.33928935898208d0*x172 + 16.0d0*x17
     &4*x175 - 4.0d0*x181) + x127 - x158
      d1E(2) = -rho_c*(x157*(x125*(x189*x212*x263 + x245*(rho_s*x213*x56
     & + x272)) + x141*(x139*x189*x266 + x195*x247*x267) + x153*(1.55214
     &407110551d0 *x186*x266 + x189*x254*x267 + x258*(x220*x255 - x225*x
     &256 + x272 )) + x184*x195*x263 - x192*x194*x23 + 0.267683468607554
     &d0*x195* x210*x266 + x203*x70*(x192*(36.0978292789602d0*x197 + 54.
     &281640924821d0*x198 + x264) + x83*(-x190 + x191)*( 16.710289428233
     &5d0*x197 + x264))) - x164*x261 + x208*x23 + x260* x28 + x262*x28)
      end subroutine


C*****************************************************************************
      pure subroutine D2ESRC_SPIN_PW92_ERF(rho_c, rho_s, mu, E, d1E, d2E
     &)
C*****************************************************************************
C   Implemented by E.R. Kjellgren.
C
C   Subroutine generated using Sympy 1.3
C   Generated: March 25, 2019
C*****************************************************************************
      implicit none
      real*8, intent(in) :: rho_c, rho_s, mu
      real*8, intent(out) :: E, d1E(9), d2E(45)
      real*8 :: x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x
     &13, x14, x15, x16, x17, x18, x19, x20, x21, x22, x23, x24, x25, x2
     &6, x27, x28, x29, x30, x31, x32, x33, x34, x35, x36, x37, x38, x39
     &, x40, x41, x42, x43, x44, x45, x46, x47, x48, x49, x50, x51, x52,
     & x53, x54, x55, x56, x57, x58, x59, x60, x61, x62, x63, x64, x65, 
     &x66, x67, x68, x69, x70, x71, x72, x73, x74, x75, x76, x77, x78, x
     &79, x80, x81, x82, x83, x84, x85, x86, x87, x88, x89, x90, x91, x9
     &2, x93, x94, x95, x96, x97, x98, x99, x100, x101, x102, x103, x104
     &, x105, x106, x107, x108, x109, x110, x111, x112, x113, x114, x115
     &, x116, x117, x118, x119, x120, x121, x122, x123, x124, x125, x126
     &, x127, x128, x129, x130, x131, x132, x133, x134, x135, x136, x137
     &, x138, x139, x140, x141, x142, x143, x144, x145, x146, x147, x148
     &, x149, x150, x151, x152, x153, x154, x155, x156, x157, x158, x159
     &, x160, x161, x162, x163, x164, x165, x166, x167, x168, x169, x170
     &, x171, x172, x173, x174, x175, x176, x177, x178, x179, x180, x181
     &, x182, x183, x184, x185, x186, x187, x188, x189, x190, x191, x192
     &, x193, x194, x195, x196, x197, x198, x199, x200, x201, x202, x203
     &, x204, x205, x206, x207, x208, x209, x210, x211, x212, x213, x214
     &, x215, x216, x217, x218, x219, x220, x221, x222, x223, x224, x225
     &, x226, x227, x228, x229, x230, x231, x232, x233, x234, x235, x236
     &, x237, x238, x239, x240, x241, x242, x243, x244, x245, x246, x247
     &, x248, x249, x250, x251, x252, x253, x254, x255, x256, x257, x258
     &, x259, x260, x261, x262, x263, x264, x265, x266, x267, x268, x269
     &, x270, x271, x272, x273, x274, x275, x276, x277, x278, x279, x280
     &, x281, x282, x283, x284, x285, x286, x287, x288, x289, x290, x291
     &, x292, x293, x294, x295, x296, x297, x298, x299, x300, x301, x302
     &, x303, x304, x305, x306, x307, x308, x309, x310, x311, x312, x313
     &, x314, x315, x316, x317, x318, x319, x320, x321, x322, x323, x324
     &, x325, x326, x327, x328, x329, x330, x331, x332, x333, x334, x335
     &, x336, x337, x338, x339, x340, x341, x342, x343, x344, x345, x346
     &, x347, x348, x349, x350, x351, x352, x353, x354, x355, x356, x357
     &, x358, x359, x360, x361, x362, x363, x364, x365, x366, x367, x368
     &, x369, x370, x371, x372, x373, x374, x375, x376, x377, x378, x379
     &, x380, x381, x382, x383, x384, x385, x386, x387, x388, x389, x390
     &, x391, x392, x393, x394, x395, x396, x397, x398, x399, x400, x401
     &, x402, x403, x404, x405, x406, x407, x408, x409, x410, x411, x412
     &, x413, x414, x415, x416, x417, x418, x419, x420, x421, x422, x423
     &, x424, x425, x426, x427, x428, x429, x430, x431, x432, x433, x434
     &, x435, x436, x437, x438, x439, x440, x441, x442, x443, x444, x445
     &, x446, x447, x448, x449, x450, x451, x452, x453, x454, x455, x456
     &, x457, x458, x459, x460, x461, x462, x463, x464, x465, x466, x467
     &, x468, x469, x470, x471, x472, x473, x474, x475, x476, x477, x478
     &, x479, x480, x481, x482, x483, x484, x485, x486, x487, x488, x489
     &, x490, x491, x492, x493, x494, x495, x496, x497, x498, x499, x500
     &, x501, x502, x503, x504, x505, x506, x507, x508, x509, x510, x511
     &, x512, x513, x514, x515, x516, x517, x518, x519, x520, x521, x522
     &, x523, x524, x525, x526, x527, x528, x529, x530, x531, x532, x533
     &, x534, x535, x536, x537, x538, x539, x540, x541, x542, x543, x544
     &, x545, x546, x547, x548, x549, x550, x551, x552, x553, x554, x555
     &, x556, x557, x558, x559, x560, x561, x562, x563, x564, x565, x566
     &, x567, x568, x569, x570, x571, x572, x573, x574, x575, x576, x577
     &, x578, x579, x580, x581, x582, x583, x584, x585, x586, x587, x588
     &, x589, x590, x591, x592, x593, x594, x595, x596, x597, x598, x599
     &, x600, x601, x602, x603, x604, x605, x606, x607, x608, x609, x610
     &, x611, x612, x613, x614, x615, x616, x617, x618, x619, x620, x621
     &, x622, x623, x624, x625, x626, x627, x628, x629, x630, x631, x632
     &, x633, x634, x635, x636, x637, x638, x639, x640, x641, x642, x643
     &, x644, x645, x646, x647, x648, x649, x650, x651, x652, x653, x654
     &, x655, x656, x657, x658, x659, x660, x661, x662, x663, x664, x665
     &, x666, x667, x668, x669, x670, x671, x672, x673, x674, x675, x676
     &, x677, x678, x679, x680, x681, x682, x683, x684, x685, x686, x687
     &, x688, x689, x690, x691, x692, x693, x694, x695, x696, x697, x698
     &, x699, x700, x701, x702, x703, x704, x705, x706, x707, x708, x709
     &, x710, x711, x712, x713, x714, x715, x716, x717, x718, x719, x720
     &, x721, x722, x723, x724, x725, x726, x727, x728, x729, x730, x731
     &, x732, x733, x734, x735, x736, x737, x738, x739, x740, x741, x742
     &, x743, x744, x745, x746, x747, x748, x749, x750, x751, x752, x753
     &, x754, x755, x756, x757, x758, x759, x760, x761, x762, x763, x764
     &, x765, x766, x767, x768, x769, x770, x771, x772, x773, x774, x775
     &, x776, x777, x778, x779, x780, x781, x782, x783, x784, x785, x786
     &, x787, x788, x789
      E = 0.0d0
      d1E(:) = 0.0d0
      d2E(:) = 0.0d0
      x0 = 2.14502939711103d0
      x1 = 1d0/x0
      x2 = rho_c**0.666666666666667d0
      x3 = 1d0/x2
      x4 = x1*x3
      x5 = 0.564189583547756d0
      x6 = rho_c**(-0.5d0)*x5
      x7 = 1.46459188756152d0
      x8 = 1d0/x7
      x9 = rho_c**(-0.333333333333333d0)*x8
      x10 = 0.826307487110758d0
      x11 = rho_c**(-0.166666666666667d0)*x10
      x12 = 7.24010193431683d0*x11 + 0.406913004517529d0*x4 + 1.41872281
     &647967d0*x6 + 3.25955091942229d0*x9
      x13 = 9.86960440108936d0
      x14 = 0.306852819440055d0
      x15 = x13/x14
      x16 = x15/x12
      x17 = log(0.5d0*x16 + 1.0d0)
      x18 = 0.194159335344114d0*x9 + 1.0d0
      x19 = 2.0d0*x18
      x20 = x17*x19
      x21 = x14/x13
      x22 = x20*x21
      x23 = 1d0/rho_c
      x24 = rho_s*x23
      x25 = x24 + 1.0d0
      x26 = -x24 + 1.0d0
      x27 = 1.92366105093154d0*x25**1.33333333333333d0 + 1.9236610509315
     &4d0*x26** 1.33333333333333d0 - 3.84732210186307d0
      x28 = rho_c**(-4)
      x29 = rho_s**4
      x30 = x28*x29
      x31 = -x30 + 1.0d0
      x32 = x27*x31
      x33 = 0.101077332976288d0*x9 + 1.0d0
      x34 = 9.87212972256927d0*x11 + 0.410025070949612d0*x4 + 0.76232752
     &193529d0*x6 + 3.29180480994506d0*x9
      x35 = 1.0d0 + 29.6085746432167d0/x34
      x36 = log(x35)
      x37 = x33*x36
      x38 = 0.0197517897025652d0*x37
      x39 = x32*x38
      x40 = 0.186690969707574d0*x9 + 1.0d0
      x41 = 13.4579137143944d0*x11 + 0.516066464547863d0*x4 + 2.91521471
     &421918d0*x6 + 5.63098414909788d0*x9
      x42 = x15/x41 + 1.0d0
      x43 = log(x42)
      x44 = x40*x43
      x45 = -x20 + x44
      x46 = x27*x30
      x47 = x21*x46
      x48 = x45*x47
      x49 = rho_c**(-2)
      x50 = rho_s**2
      x51 = x49*x50
      x52 = -x51 + 1.0d0
      x53 = 0.0837224526465878d0
      x54 = exp(-0.6836107611867116d0*x9)
      x55 = 4.60115111447049d0
      x56 = 1d0/x55
      x57 = rho_c**(-1.33333333333333d0)
      x58 = x56*x57
      x59 = rho_c**1.0d0
      x60 = 1d0/x59
      x61 = 3.14159265358979d0
      x62 = 1d0/x61
      x63 = 0.009578475d0*x62
      x64 = 0.0676322201645715d0*x4 + 0.00126674656487366d0*x58 - x60*x6
     &3 + 0.0188171922990732d0*x9 + 1.0d0
      x65 = x54*x64
      x66 = x53*x65
      x67 = 1.41421356237310d0
      x68 = mu**5*x67
      x69 = x66*x68
      x70 = rho_c**(-1.66666666666667d0)
      x71 = 0.0892278228691848d0*x70
      x72 = x25**0.666666666666667d0
      x73 = x26**0.666666666666667d0
      x74 = x72 + x73
      x75 = x74**3
      x76 = mu**2
      x77 = x74**2
      x78 = 1d0/x77
      x79 = x76*x78
      x80 = x79*x9
      x81 = 1d0/x74
      x82 = mu*x81
      x83 = 11.1447260721495d0*x11*x82 + 1.0d0
      x84 = 12.5327170711751d0*x80 + x83
      x85 = 1d0/x84
      x86 = mu**3
      x87 = 1d0/x75
      x88 = 27.1408204624105d0*x6*x86*x87 + 27.0733719592201d0*x80 + x83
     &
      x89 = x85*x88
      x90 = x21*log(x89)
      x91 = 0.179587122125167d0
      x92 = x65*x91
      x93 = 0.175432441092201d0*x60
      x94 = x52*x93
      x95 = 0.825481812223657d0*x4 - 4.49737346725955d0*x9
      x96 = exp(-0.28165369188898165d0*x9)
      x97 = x61*x96
      x98 = x52*x97
      x99 = 12.0d0*x98
      x100 = x95*x99
      x101 = 2.28942848510666d0
      x102 = x101*x55
      x103 = 0.817704266774463d0*x102
      x104 = x25**2
      x105 = 1d0/x25
      x106 = x105**0.333333333333333d0
      x107 = x106*x9
      x108 = 0.0259335011650457d0*x107
      x109 = -x108 + 1.0d0
      x110 = x105**0.666666666666667d0
      x111 = 0.0524148278841779d0*x4
      x112 = 0.494402081358784d0*x107 + x110*x111 + 1.0d0
      x113 = 1d0/x112
      x114 = 1d0/x110
      x115 = x113*x114
      x116 = x109*x115
      x117 = x104*x116
      x118 = x103*x117
      x119 = x26**2
      x120 = 1d0/x26
      x121 = x120**0.333333333333333d0
      x122 = x121*x9
      x123 = 0.0259335011650457d0*x122
      x124 = -x123 + 1.0d0
      x125 = x120**0.666666666666667d0
      x126 = x111*x125 + 0.494402081358784d0*x122 + 1.0d0
      x127 = 1d0/x126
      x128 = 1d0/x125
      x129 = x127*x128
      x130 = x124*x129
      x131 = x119*x130
      x132 = x103*x131
      x133 = x53*(x100*x59 + x118*x2 + x132*x2)
      x134 = 0.00111534778586481d0*x70
      x135 = x67*x86
      x136 = 4.0d0*x135
      x137 = -x22 + x39
      x138 = x137 - x48
      x139 = 0.0472353356922751d0*mu**8
      x140 = x138*x139
      x141 = rho_c**(-2.66666666666667d0)
      x142 = 0.267683468607554d0*x141
      x143 = 0.148394183600699d0
      x144 = x65 - 1.0d0
      x145 = x143*x144
      x146 = x145*x52
      x147 = rho_c**(-2.0d0)
      x148 = 0.101321183642338d0
      x149 = 0.526297323276602d0*x148
      x150 = x147*x149
      x151 = mu**6
      x152 = 4.0d0*x151
      x153 = x144*x52
      x154 = 0.13157433081915d0*x153
      x155 = x60*x62
      x156 = x102*(x25**2.66666666666667d0 + x26**2.66666666666667d0)
      x157 = 0.817704266774463d0*x156
      x158 = exp(-0.49698248213959023d0*x9)
      x159 = x0*x158
      x160 = 0.558025705063192d0*x4 - 0.352521395009435d0*x9
      x161 = x160*x52
      x162 = x159*x161
      x163 = 10.9027235569928d0*x162
      x164 = x118 + x132 - x157 + x163
      x165 = 0.00139418473233101d0*x143
      x166 = 1.55214407110551d0*x58
      x167 = mu**4
      x168 = 4.0d0*x167
      x169 = -x136*(x133*x134 + x92*x94) + x140*x142 + x152*(x138*x150 -
     & 0.0167302167879722d0*x146*x70) + x168*(-rho_c**(-1.0d0)*x164*x165
     & + x138*x166 - x154*x155) - x52*x69*x71 - x75*x90
      x170 = 0.508616435555896d0*x4*x76 + 1.0d0
      x171 = x170**(-4)
      x172 = 0.25d0*x171
      x173 = x169*x172
      x174 = 0.0691804894611514d0*x14
      x175 = x17*x174
      x176 = x175*x57
      x177 = x57*x8
      x178 = x32*x36
      x179 = x177*x178
      x180 = rho_c**(-5)
      x181 = x180*x29
      x182 = x27*x37
      x183 = x181*x182
      x184 = x26**0.333333333333333d0
      x185 = x25**0.333333333333333d0
      x186 = x184 - x185
      x187 = x31*x37
      x188 = rho_s*x49
      x189 = x186*x187*x188
      x190 = x1*x70
      x191 = rho_c**(-1.5d0)*x5
      x192 = rho_c**(-1.16666666666667d0)*x10
      x193 = 1.08651697314076d0*x177 + 0.27127533634502d0*x190 + 0.70936
     &1408239834d0* x191 + 1.20668365571947d0*x192
      x194 = x16 + 2.0d0
      x195 = 1d0/x194
      x196 = x195/x12**2
      x197 = x193*x196
      x198 = x18*x197
      x199 = 1.09726826998169d0*x177 + 0.273350047299742d0*x190 + 0.3811
     &63760967645d0 *x191 + 1.64535495376155d0*x192
      x200 = 1d0/x35
      x201 = x34**(-2)
      x202 = x32*x33
      x203 = x200*x201*x202
      x204 = x199*x203
      x205 = x21*x45
      x206 = rho_s**5
      x207 = rho_c**(-6)
      x208 = x206*x207
      x209 = x186*x208
      x210 = x205*x209
      x211 = 16.0d0*x205
      x212 = x181*x27
      x213 = 1.87699471636596d0*x177 + 0.344044309698576d0*x190 + 1.4576
     &0735710959d0* x191 + 2.24298561906574d0*x192
      x214 = 1d0/x42
      x215 = x41**(-2)
      x216 = x214*x215*x40
      x217 = x213*x216
      x218 = 0.062230323235858d0*x57
      x219 = x174*x43
      x220 = x218*x219
      x221 = x19*x197
      x222 = 0.129439556896076d0*x57
      x223 = x175*x222
      x224 = -x221 + x223
      x225 = x217 - x220 + x224
      x226 = x225*x46
      x227 = x170**(-5)
      x228 = x169*x227*x76
      x229 = x190*x228
      x230 = rho_c**(-4.66666666666667d0)
      x231 = x230*x50
      x232 = 0.17845564573837d0*x69
      x233 = x52*x54
      x234 = rho_c**(-2.33333333333333d0)
      x235 = x234*x56
      x236 = 0.00168899541983155d0*x235
      x237 = x147*x63
      x238 = 0.045088146776381d0*x190
      x239 = 0.00627239743302441d0*x177
      x240 = x236 - x237 + x238 + x239
      x241 = x240*x53*x68
      x242 = x65*x68
      x243 = x242*x52
      x244 = rho_c**(-3.0d0)
      x245 = 0.0571643564037363d0
      x246 = x244*x245
      x247 = x141*x53
      x248 = x25**(-0.333333333333333d0)
      x249 = x26**(-0.333333333333333d0)
      x250 = x248 - x249
      x251 = rho_s*x250
      x252 = 2.0d0*x77*x90
      x253 = x252*x49
      x254 = rho_c**(-3.66666666666667d0)
      x255 = 1.85745434535825d0*x192
      x256 = x177*x82
      x257 = 9.02445731974005d0*x256
      x258 = x191*x79
      x259 = 13.5704102312052d0*x258
      x260 = rho_c**(-2.16666666666667d0)*x10
      x261 = x260*x81
      x262 = x251*x261
      x263 = 7.42981738143299d0*x262
      x264 = rho_c**(-2.33333333333333d0)
      x265 = rho_s*x264
      x266 = mu*x78
      x267 = x266*x8
      x268 = x250*x265*x267
      x269 = 36.0978292789602d0*x268
      x270 = rho_c**(-2.5d0)*x5
      x271 = x270*x76
      x272 = x251*x87
      x273 = 54.281640924821d0*x271*x272
      x274 = x255 - x263
      x275 = 4.17757235705837d0*x256 - 16.7102894282335d0*x268 + x274
      x276 = x21*(-x255 - x257 - x259 + x263 + x269 + x273 + x275*x89)
      x277 = 1d0/x88
      x278 = mu*x277
      x279 = x278*x77
      x280 = 4.0d0*x27
      x281 = x181*x280
      x282 = 0.0506609980493537d0*x187
      x283 = x186*x282
      x284 = -0.000665486074881265d0*x179 + 0.0790071588102608d0*x183 + 
     &x188*x283 + 0.58482233974552d0*x204 + x224
      x285 = x205*x281 - 2.56488140124205d0*x210 - x226 + x284
      x286 = x139*x285
      x287 = 0.350864882184401d0*x92
      x288 = rho_c**(-4.0d0)
      x289 = x288*x50
      x290 = x54*x91
      x291 = x240*x290
      x292 = 0.122619224952946d0
      x293 = x292*x65
      x294 = 0.00185891297644135d0*x141
      x295 = 0.550321208149104d0*x190
      x296 = 1.49912448908652d0*x177
      x297 = -x295 + x296
      x298 = x59*x99
      x299 = x297*x298
      x300 = x95*x97
      x301 = 24.0d0*x300
      x302 = x301*x50
      x303 = x147*x302
      x304 = x52*x95*x96
      x305 = 2.14502939711103d0
      x306 = rho_c**(-0.333333333333333d0)*x305
      x307 = 1.12661476755593d0*x304*x306
      x308 = 0.00706864485168613d0*x2
      x309 = x105*x265
      x310 = -x309 + x57
      x311 = x101*x61
      x312 = 1d0/x106
      x313 = x113*x312
      x314 = x311*x313
      x315 = x104*x310*x314
      x316 = x308*x315
      x317 = 1d0/x121
      x318 = x127*x317
      x319 = x311*x318
      x320 = x119*x319
      x321 = x120*x265
      x322 = x321 + x57
      x323 = x308*x322
      x324 = x320*x323
      x325 = rho_c**(-0.333333333333333d0)
      x326 = x102*x325
      x327 = 0.545136177849642d0*x326
      x328 = x117*x327
      x329 = x131*x327
      x330 = rho_c**(-1.33333333333333d0)
      x331 = x102*x330
      x332 = rho_s*x331
      x333 = 2.18054471139857d0*x332
      x334 = x130*x26
      x335 = x333*x334
      x336 = x112**(-2)
      x337 = x114*x336
      x338 = x109*x337
      x339 = x103*x2
      x340 = x338*x339
      x341 = x1*x141
      x342 = 0.034943218589452d0*x341
      x343 = rho_s*x105
      x344 = x110*x343
      x345 = x342*x344
      x346 = 0.164800693786261d0*x8
      x347 = x106*x309
      x348 = x346*x347
      x349 = 0.034943218589452d0*x190
      x350 = x110*x349
      x351 = 0.164800693786261d0*x177
      x352 = x106*x351
      x353 = x350 + x352
      x354 = -x345 - x348 + x353
      x355 = x104*x354
      x356 = x340*x355
      x357 = x126**(-2)
      x358 = x128*x357
      x359 = x124*x358
      x360 = x339*x359
      x361 = rho_s*x120
      x362 = x125*x361
      x363 = x121*x321
      x364 = x121*x351 + x125*x349
      x365 = x342*x362 + x346*x363 + x364
      x366 = x119*x365
      x367 = x360*x366
      x368 = x116*x25
      x369 = x333*x368
      x370 = x316 + x324 + x328 + x329 + x335 + x356 + x367 - x369
      x371 = x134*x53
      x372 = 0.0334604335759443d0*x145
      x373 = x141*x143
      x374 = 0.000844497709915773d0*x235
      x375 = x147*x62
      x376 = 0.0047892375d0*x375
      x377 = 0.0225440733881905d0*x190
      x378 = 0.0031361987165122d0*x177
      x379 = x54*(0.113935126864452d0*x177*x64 - x374 + x376 - x377 - x3
     &78)
      x380 = x379*x52
      x381 = 0.0334604335759443d0*x143*x70
      x382 = x20 - x44
      x383 = x137 + x382*x47
      x384 = x148*x244
      x385 = -x217 + x220 + x221 - x223
      x386 = x21*x382
      x387 = 2.56488140124205d0*x386
      x388 = x209*x387 - x281*x386 + x284 + x385*x46
      x389 = 0.263148661638301d0*x144
      x390 = x50*x62
      x391 = x288*x390
      x392 = 0.263148661638301d0*x155
      x393 = 0.00232364122055169d0*x143
      x394 = 0.545136177849642d0*x156*x325
      x395 = 7.26848237132856d0*x162*x325
      x396 = 0.117507131669812d0*x177 - 0.372017136708795d0*x190
      x397 = x159*x396
      x398 = x397*x52
      x399 = 10.9027235569928d0*x2
      x400 = x398*x399
      x401 = x159*x160
      x402 = 21.8054471139857d0*x401
      x403 = x264*x50
      x404 = x402*x403
      x405 = x158*x161
      x406 = x3*x7
      x407 = 1.80615420514536d0*x405*x406
      x408 = x26**1.66666666666667d0
      x409 = x25**1.66666666666667d0
      x410 = x102*(x408 - x409)
      x411 = rho_s*x330
      x412 = 2.18054471139857d0*x410*x411
      x413 = -x316 - x324 - x328 - x329 - x335 - x356 - x367 + x369
      x414 = x165*x70
      x415 = x136*(x133*x294 + 0.175432441092201d0*x147*x52*x92 - 0.0399
     &758348639607d0*x234*x293*x52 - x287*x289 + x291*x94 - x371* (x100 
     &+ x299 + x303 + x307 + x370)) - 0.713822582953479d0*x140* x254 + x
     &142*x286 - x152*(-x150*x388 - 0.0278836946466203d0*x153* x373 + x2
     &31*x372 + x380*x381 + 1.0525946465532d0*x383*x384) - x168*(x147*x3
     &93*(-x118 - x132 + x157 - x163) - x154*x375 - x166* x388 + 2.06952
     &542814068d0*x235*x383 + x380*x392 + x389*x391 - x414*(x394 - x395 
     &- x400 - x404 - x407 + x412 + x413)) - x231* x232 + x233*x241*x71 
     &- 0.0203323666368788d0*x243*x246 + 0.148713038115308d0*x243*x247 +
     & x251*x253 - x276*x279
      x416 = x171*x415
      x417 = rho_s**3
      x418 = x28*x417
      x419 = 0.0790071588102608d0*x182
      x420 = x418*x419
      x421 = x23*x283
      x422 = 2.56488140124205d0*x181
      x423 = x186*x205
      x424 = x422*x423
      x425 = x280*x418
      x426 = x205*x425
      x427 = rho_s*x254
      x428 = x23*x250
      x429 = x278*x74
      x430 = 7.42981738143299d0*x192
      x431 = 36.0978292789602d0*x256 + 54.281640924821d0*x258 + x430
      x432 = x250*x431
      x433 = -x248 + x249
      x434 = 16.7102894282335d0*x256 + x430
      x435 = x433*x434
      x436 = x21*(x432 + x435*x89)
      x437 = rho_c**(-3)
      x438 = x417*x437
      x439 = -x419*x438
      x440 = x280*x438
      x441 = -x205*x440 - x283 + 2.56488140124205d0*x30*x423 + x439
      x442 = x139*x441
      x443 = x254*x372
      x444 = x244*x287
      x445 = x301*x60
      x446 = x108 - 1.0d0
      x447 = x115*x446
      x448 = 2.18054471139857d0*x326
      x449 = x25*x448
      x450 = x26*x448
      x451 = x123 - 1.0d0
      x452 = x129*x451
      x453 = x25*x339
      x454 = x337*x446
      x455 = x353*x454
      x456 = x26*x339
      x457 = x358*x451
      x458 = x364*x457
      x459 = x26*x318
      x460 = x3*x311
      x461 = 0.00706864485168613d0*x460
      x462 = x25*x313
      x463 = x459*x461 - x461*x462
      x464 = x447*x449 - x450*x452 + x453*x455 - x456*x458 + x463
      x465 = rho_s*x445 + x464
      x466 = x244*x62
      x467 = x389*x466
      x468 = x330*x402
      x469 = rho_s*x468
      x470 = 2.18054471139857d0*x325
      x471 = x136*(rho_s*x444 + x371*x465) + x152*(rho_s*x443 + x149*x24
     &4*x441) + x168*(rho_s*x467 + 1.55214407110551d0*x235*x441 + x414*(
     &-x410* x470 + x464 + x469)) + x232*x427 - x252*x428 + 0.2676834686
     &07554d0*x254*x442 + x429*x436
      x472 = x172*x471
      x473 = x175*x234
      x474 = x27*x29
      x475 = x36*x8
      x476 = rho_c**(-6.33333333333333d0)*x474*x475
      x477 = x206/rho_c**7
      x478 = x186*x37
      x479 = 4.86345581273796d0*x478
      x480 = x186*x31
      x481 = rho_c**(-3.33333333333333d0)
      x482 = rho_s*x481
      x483 = x475*x480*x482
      x484 = x234*x8
      x485 = x178*x484
      x486 = x207*x29
      x487 = x182*x486
      x488 = 0.237021476430782d0*x187
      x489 = x25**(-0.666666666666667d0)
      x490 = 0.854960467080683d0*x24
      x491 = x489*x490
      x492 = x26**(-0.666666666666667d0)
      x493 = x490*x492
      x494 = x491 + x493
      x495 = -5.1297628024841d0*x184 + 5.1297628024841d0*x185 + x494
      x496 = rho_s*x437
      x497 = x495*x496
      x498 = 0.682784063255296d0
      x499 = x197*x498
      x500 = x177*x197
      x501 = rho_c**(-2.16666666666667d0)*x10
      x502 = x196*(1.06404211235975d0*x270 + 0.452125560575033d0*x341 + 
     &1.44868929752102d0*x484 + 1.40779759833938d0*x501)
      x503 = 24.0d0*x18
      x504 = x15*x193**2/(x12**4*x194**2)
      x505 = x199*x200
      x506 = x201*x505
      x507 = x177*x32*x506
      x508 = x33*x506
      x509 = x212*x508
      x510 = x203*(0.571745641451468d0*x270 + 0.455583412166236d0*x341 +
     & 1.46302435997558d0*x484 + 1.91958077938847d0*x501)
      x511 = x188*x480*x508
      x512 = x199**2*x202/(x34**4*x35**2)
      x513 = x193*x195*(2.17303394628153d0*x177 + 0.542550672690039d0*x1
     &90 + 1.41872281647967d0*x191 + 2.41336731143894d0*x192)/x12**3
      x514 = x205*x477
      x515 = x186*x514
      x516 = x205*x207*x474
      x517 = x202*x505*(2.19453653996337d0*x177 + 0.546700094599483d0*x1
     &90 + 0.76232752193529d0*x191 + 3.29070990752309d0*x192)/x34**3
      x518 = x495*x514
      x519 = x186*x225
      x520 = x208*x519
      x521 = x225*x27
      x522 = x181*x521
      x523 = 0.172586075861435d0*x473
      x524 = x213*x214
      x525 = x215*x524
      x526 = x222*x499
      x527 = 0.129439556896076d0*x500
      x528 = x19*x502
      x529 = x19*x504
      x530 = x19*x513
      x531 = x46*(x15*x213**2*x40/(x41**4*x42**2) + 0.062230323235858d0*
     &x177*x525 + x216*(2.18641103566438d0*x270 + 0.573407182830959d0*x3
     &41 + 2.50265962182128d0*x484 + 2.6168165555767d0*x501) + x218*x498
     &* x525 - 0.0829737643144773d0*x219*x234 - x40*x524*( 3.75398943273
     &192d0*x177 + 0.688088619397151d0*x190 + 2.91521471421918d0*x191 + 
     &4.48597123813148d0*x192)/x41**3 + x523 - x526 - x527 - x528 - x529
     & + x530)
      x532 = x481*x56
      x533 = x190*x227*x76
      x534 = x241*x54
      x535 = 0.121994199821273d0*x242*x245
      x536 = x50*x535
      x537 = rho_c**(-5.66666666666667d0)*x50
      x538 = x233*x68
      x539 = x240*x538
      x540 = 0.0751469112939684d0*x341 - 0.01915695d0*x466 + 0.008363196
     &57736588d0* x484 + 0.00394098931294027d0*x532
      x541 = rho_c**(-4.33333333333333d0)
      x542 = x254*x69
      x543 = x28*x50
      x544 = x250**2
      x545 = 8.0d0*x544*x74*x90
      x546 = x25**(-1.33333333333333d0)
      x547 = x24*x546
      x548 = x26**(-1.33333333333333d0)
      x549 = x24*x548
      x550 = x547 + x549
      x551 = rho_s*(-6.0d0*x248 + 6.0d0*x249 + x550)
      x552 = x251*x49
      x553 = x276*x429
      x554 = 3.0d0*x76
      x555 = x276*x554
      x556 = x277*x85
      x557 = x555*x556
      x558 = x257 + x259 - x269 - x273 + x274
      x559 = x88**(-2)
      x560 = x555*x559
      x561 = x484*x82
      x562 = 5.57009647607783d0*x561
      x563 = x544*x87
      x564 = mu*x50*x541*x563*x8
      x565 = rho_c**(-3.33333333333333d0)
      x566 = x267*x551*x565
      x567 = x251*x267*x481
      x568 = x10*x81
      x569 = 2.476605793811d0*x568
      x570 = 9.90642317524398d0*x544
      x571 = rho_c**(-3.16666666666667d0)
      x572 = x251*x571
      x573 = rho_c**(-4.16666666666667d0)*x10*x50*x570*x78 + rho_c**( -3
     &.16666666666667d0)*x551*x569 + 2.16703006958462d0*x501 - x569* x57
     &2
      x574 = 12.0326097596534d0*x561
      x575 = x270*x79
      x576 = x275*x82
      x577 = 2.0d0*x85
      x578 = x88/x84**2
      x579 = x576*x578
      x580 = x5*x76
      x581 = rho_c**(-3.5d0)*x580
      x582 = 3.0d0*x21
      x583 = 0.40528798439483d0*x478
      x584 = x31*x38
      x585 = -0.00532388859905012d0*x476 + x477*x583 - 0.003413785712497
     &06d0*x483 + 0.000887314766508353d0*x485 - 0.395035794051304d0*x487
     & + x497* x584 - 0.0394081882442864d0*x507 + 4.67857871796416d0*x50
     &9 - 0.58482233974552d0*x510 + 2.99999988448829d0*x511 - 17.3157558
     &993759d0*x512 + 20.5190512099364d0*x515 - 20.0d0*x516 + 0.58482233
     &974552d0*x517 - x518 - 5.1297628024841d0*x520 + 8.0d0* x522 - x523
     & + x526 + x527 + x528 + x529 - x530 + x531
      x586 = 0.803050405822664d0*x139
      x587 = x141*x586
      x588 = -x236 + x237 - x238 - x239
      x589 = x290*x588
      x590 = rho_c**(-5.33333333333333d0)
      x591 = rho_c**(-5.0d0)
      x592 = x51 - 1.0d0
      x593 = x305*x95
      x594 = 48.0d0*x50
      x595 = 2.25322953511185d0*x96
      x596 = 3.63424118566428d0*x116
      x597 = x50*x565
      x598 = x102*x597
      x599 = x596*x598
      x600 = x311*x462
      x601 = 0.0376994392089927d0*x411
      x602 = x310*x600*x601
      x603 = 0.00942485980224818d0*x325
      x604 = x315*x603
      x605 = 3.63424118566428d0*x130
      x606 = x598*x605
      x607 = x102*x265
      x608 = 1.45369647426571d0*x607
      x609 = x368*x608
      x610 = 0.181712059283214d0*x331
      x611 = x117*x610
      x612 = x311*x459
      x613 = x322*x612
      x614 = x601*x613
      x615 = x320*x322*x603
      x616 = x334*x608
      x617 = x131*x610
      x618 = 0.0115260005177981d0*x234
      x619 = 0.0201705009061467d0*x565
      x620 = 0.00288150012944953d0*x481
      x621 = x50/x104
      x622 = rho_c**(-4.33333333333333d0)
      x623 = 0.0115260005177981d0*x622
      x624 = 0.817704266774463d0*x2
      x625 = x104*x314*x624*(x343*x619 + x343*x620 - x618 - x621*x623)
      x626 = x50/x119
      x627 = x320*x624*(x361*x619 + x361*x620 + x618 + x623*x626)
      x628 = x312*x336
      x629 = x354*x628
      x630 = x310*x311
      x631 = 0.0141372897033723d0*x2
      x632 = x104*x629*x630*x631
      x633 = x338*x354
      x634 = 4.36108942279713d0*x332
      x635 = x25*x633*x634
      x636 = 1.09027235569928d0*x326
      x637 = x338*x355*x636
      x638 = x317*x357
      x639 = x365*x638
      x640 = x119*x311*x322*x631*x639
      x641 = x26*x359*x365*x634
      x642 = x359*x366*x636
      x643 = 0.0582386976490866d0*x341
      x644 = 0.219734258381682d0*x484
      x645 = x106*x644
      x646 = x1*x254
      x647 = 0.116477395298173d0*x646
      x648 = x106*x8
      x649 = x343*x648
      x650 = 0.384534952167943d0*x565
      x651 = 0.0549335645954204d0*x481
      x652 = 0.0582386976490866d0*x1*x230
      x653 = 0.219734258381682d0*x622
      x654 = x104*x340*(-x110*x621*x652 - x110*x643 + x344*x647 - x621*x
     &648*x653 - x645 + x649*x650 + x649*x651)
      x655 = x112**(-3)
      x656 = 0.0698864371789039d0*x341
      x657 = 0.329601387572523d0*x8
      x658 = 0.0698864371789039d0*x190
      x659 = 0.329601387572523d0*x177
      x660 = x106*x659 + x110*x658
      x661 = x109*x339
      x662 = x114*x355*x655*x661*(-x344*x656 - x347*x657 + x660)
      x663 = x121*x644
      x664 = x121*x8
      x665 = x361*x664
      x666 = x119*x360*(x125*x626*x652 + x125*x643 + x362*x647 + x626*x6
     &53*x664 + x650*x665 + x651*x665 + x663)
      x667 = x126**(-3)
      x668 = x121*x659 + x125*x658
      x669 = x124*x339
      x670 = x128*x366*x667*x669*(x362*x656 + x363*x657 + x668)
      x671 = 12.0d0*x135
      x672 = x143*x379
      x673 = x233*(-0.455740507457808d0*x177*(x374 - x376 + x377 + x378)
     & - x244*x63 + 0.0259624262672375d0*x341*x64 + 0.0375734556469842d0
     &*x341 - 0.151913502485936d0*x484*x64 + 0.00418159828868294d0*x484 
     &+ 0.00197049465647014d0*x532)
      x674 = x148*x288
      x675 = 12.0d0*x151
      x676 = 3.61230841029072d0*x158
      x677 = x160*x676*x7
      x678 = x50*x677
      x679 = 5.0d0*x24
      x680 = x679*x72
      x681 = x679*x73
      x682 = 12.0d0*x167
      x683 = 0.0833333333333333d0*rho_c
      x684 = x27*x417
      x685 = x475*x590*x684
      x686 = x36*x484
      x687 = x180*x417
      x688 = x182*x687
      x689 = 2.56488140124205d0*x185
      x690 = 2.56488140124205d0*x184
      x691 = x494 + x689 - x690
      x692 = x27*x418*x508
      x693 = x23*x508
      x694 = x205*x486
      x695 = x186*x694
      x696 = x180*x684
      x697 = x691*x694
      x698 = x418*x521
      x699 = rho_s*x230
      x700 = 3.0d0*x248
      x701 = 3.0d0*x249
      x702 = x550 - x700 + x701
      x703 = x278*x436
      x704 = x431*x433
      x705 = x250*x434
      x706 = 0.00170689285624853d0*x686
      x707 = x49*x584
      x708 = 1.49999994224414d0*x693
      x709 = 0.00266194429952506d0*x685 + 0.316028635241043d0*x688 - 2.3
     &3928935898208d0*x692
      x710 = 1.2383028969055d0*x501
      x711 = 72.1956585579204d0*x544
      x712 = 144.751042466189d0*x563
      x713 = 2.476605793811d0*x260
      x714 = x264*x8
      x715 = x714*x82
      x716 = 12.0326097596534d0*x715
      x717 = 18.0938803082737d0*x575
      x718 = 14.859634762866d0*x192 + 33.420578856467d0*x256
      x719 = -x547 - x549 + x700 - x701
      x720 = 5.57009647607783d0*x715
      x721 = x429*x582
      x722 = -x184 + x185
      x723 = x282*x722 + x30*x387*x722 + x386*x440 + x439
      x724 = -x491 - x493 - x689 + x690
      x725 = 0.202643992197415d0*x486
      x726 = x31*x722
      x727 = x386*x486
      x728 = 10.2595256049682d0*x727
      x729 = x186*x728 + x37*x722*x725 + x385*x422*x722 + x385*x425 - 16
     &.0d0*x386* x696 - x478*x725 - x706*x726 + x707*x724 + x708*x726 + 
     &x709 - x722*x728 + x724*x727
      x730 = rho_s*x288
      x731 = 0.726848237132856d0*x331
      x732 = x26*x731
      x733 = x309 - x57
      x734 = 0.0188497196044964d0*x325
      x735 = x600*x734
      x736 = x345 + x348 - x350 - x352
      x737 = x25*x461
      x738 = 0.0115260005177981d0*x565
      x739 = x343*x738
      x740 = 0.00288150012944953d0*x234
      x741 = 0.00864450038834858d0*x264
      x742 = x740 + x741
      x743 = x600*x624
      x744 = 3.63424118566428d0*x607
      x745 = x353*x628
      x746 = x25*x308*x745
      x747 = x25*x327
      x748 = 0.0582386976490866d0*x341
      x749 = x125*x748
      x750 = x264*x346
      x751 = 0.0549335645954204d0*x484
      x752 = 0.0582386976490866d0*x646
      x753 = 0.219734258381682d0*x565
      x754 = x121*x750 + x121*x751 + x362*x752 + x665*x753 + x749
      x755 = x26*x327
      x756 = x365*x450
      x757 = x110*x748
      x758 = x106*x750
      x759 = x106*x751
      x760 = x344*x752
      x761 = x649*x753
      x762 = x114*x655*x660
      x763 = x453*x762
      x764 = x128*x667*x668
      x765 = x365*x456*x764
      x766 = x311*x70
      x767 = 0.00471242990112409d0*x766
      x768 = 0.0188497196044964d0*rho_c**(-2.66666666666667d0)*rho_s
      x769 = x364*x638
      x770 = x26*x311*x323*x769 + x26*x461*x639 + x314*x768 + x319*x768 
     &+ x459*x767 - x462*x767 - x612*x624*(x361*x738 + x742) + x613*x734
     &
      x771 = x102*(-x408 + x409)
      x772 = x359*x364
      x773 = x338*x353
      x774 = x489 + x492
      x775 = 0.0168869993497846d0*x187*x774
      x776 = 20.5190512099364d0*x423
      x777 = 12.0d0*x205*x27
      x778 = x546 + x548
      x779 = x436*x554*x81
      x780 = -0.237021476430782d0*x182*x51 - 0.854960467080683d0*x205*x3
     &0*x774 + x438 *x583 + x438*x776 - x51*x777 + x775
      x781 = x713*x778
      x782 = 9.90642317524398d0*x261
      x783 = x266*x714
      x784 = x433**2
      x785 = 0.0282745794067445d0*x766
      x786 = 0.0141372897033723d0*x460
      x787 = 3.63424118566428d0*x331
      x788 = 4.36108942279713d0*x326
      x789 = x116*x787 + x130*x787 + x313*x785 + x318*x785 - x340*(x645 
     &+ x757) + x353*x661*x762 - x360*(x663 + x749) + x364*x669*x764 + x
     &745*x786 + x769*x786 + x772*x788 + x773*x788
      E = -rho_c*(x173 + x22 - x39 + x48)
      d1E(1) = 0.25d0*rho_c*(0.517758227584304d0*x176 - 0.00266194429952
     &506d0*x179 + 0.316028635241043d0*x183 + 0.202643992197415d0*x189 -
     & 8.0d0*x198 + 2.33928935898208d0*x204 - 10.2595256049682d0*x210 + 
     &x211*x212 - 4.0d0*x226 - 1.35631049481572d0*x229 - x416) + x138 - 
     &x173
      d1E(2) = -rho_c*(x420 + x421 - x424 + x426 + x472)
      d2E(1) = 0.258879113792152d0*x176 - 0.00133097214976253d0*x179 + 0
     &.158014317620522d0*x183 + 0.101321996098707d0*x189 - 4.0d0*x198 + 
     &1.16964467949104d0*x204 + 8.0d0*x205*x212 - 5.1297628024841d0* x21
     &0 - 2.0d0*x226 - 0.678155247407862d0*x229 - 0.5d0*x416 + x683* (-6
     &.89841809380227d0*x167*x169*x532/x170**6 - x171*(-rho_c**( -6.0d0)
     &*x536 - rho_c**(-6.0d0)*x536 + 7.85204841248826d0*x140* x230 + 1.0
     &7073387443022d0*x231*x534 + 0.284653132916304d0*x243* x245*x288 - 
     &0.000542507213303895d0*x243*x541 + 0.121994199821273d0*x246*x539 -
     & 0.892278228691848d0*x247*x539 + x252*x437*x551 - 4.28293549772087
     &d0*x254*x286 + x275*x557*x74 - x279*x582*(144.751042466189d0*rho_c
     &**(-4.5d0)*x50*x544*x580/x74** 4 - 54.281640924821d0*x272*x581 + 1
     &8.0938803082737d0*x551*x581* x87 - x558*x576*x577 + 72.19565855792
     &04d0*x564 + 12.0326097596534d0*x566 - 24.0652195193068d0*x567 + x5
     &73 + x574 + 20.3556153468079d0*x575 + x579*(3.71490869071649d0*x19
     &2 + 8.35514471411675d0*x256 - 14.859634762866d0*x262 - 33.42057885
     &6467d0*x268) - x89*(x562 + 33.420578856467d0*x564 + 5.570096476077
     &83d0*x566 - 11.1401929521557d0*x567 + x573)) - 1.18970430492246d0*
     &x52*x542 - 0.267683468607554d0*x53*x538*x540* x70 + 3.390657269029
     &02d0*x537*x69 - x543*x545 + 12.0d0*x552*x553 - x558*x560*x74 + x58
     &5*x587 - x671*(-0.00910930363347549d0*rho_c **(-3.66666666666667d0
     &)*x592*x66 + 0.00495710127051027d0*x133* x254 + 0.350864882184401d
     &0*x147*x589*x592 - 0.0799516697279215d0* x234*x292*x54*x588*x592 +
     & 0.0037178259528827d0*x247*(-x100 - x299 - x303 - x307 + x413) + 0
     &.701729764368802d0*x289*x589 - x290*x540 *x592*x93 + 0.13325278287
     &9869d0*x293*x481*x592 + 0.159903339455843d0*x293*x50*x590 - x371*(
     &-x147*x297*x594*x97 + x244*x300*x594 - x297*x306*x52*x595 - 24.0d0
     &*x297*x98 - x298*( 0.917202013581841d0*x341 - 1.99883265211535d0*x
     &484) - x302*x437 - 0.751076511703951d0*x304*x305*x57 - 0.154912426
     &780983d0*x304*x70 - 4.50645907022371d0*x481*x50*x593*x96 - x599 + 
     &x602 - x604 - x606 - x609 + x611 - x614 - x615 + x616 + x617 - x62
     &5 + x627 - x632 + x635 - x637 - x640 - x641 - x642 - x654 - x662 +
     & x666 - x670) - x444*x592 - 1.75432441092201d0*x50*x591*x92) + x67
     &5*( 3.15778393965961d0*x138*x674 + 0.211916079314314d0*x145*x537 -
     & 0.074356519057654d0*x146*x254 + x150*x585 - 0.133841734303777d0* 
     &x231*x672 - 2.10518929310641d0*x285*x384 + 0.111534778586481d0* x3
     &73*x380 - x381*x673) + x682*(4.82889266566159d0*x138*x532 - 0.0061
     &9637658813783d0*x143*x164*x244 + 1.3157433081915d0*x144* x390*x591
     & + x166*x585 - 4.13905085628136d0*x235*x285 + 0.00464728244110338d
     &0*x373*(x370 - x394 + x395 + x400 + x404 + x407 - x412) + 0.526297
     &323276602d0*x375*x380 - 1.0525946465532d0* x379*x391 - x392*x673 -
     & x414*(rho_c**(-3.66666666666667d0)*x678 + 0.299209d0*x147*x405 + 
     &0.181712059283214d0*x156*x330 + x159*x399* x52*(0.620028561181324d
     &0*x341 - 0.156676175559749d0*x484) - 2.42282745710952d0*x162*x330 
     &+ x254*x678 - 2.90739294853142d0* x265*x410 + 14.5369647426571d0*x
     &325*x398 + x396*x406*x52*x676 + 43.6108942279714d0*x397*x403 - 36.
     &3424118566428d0*x401*x597 + x599 - x602 + x604 + x606 - 0.72684823
     &7132856d0*x607*(-6.0d0*x408 + 6.0d0*x409 + x680 + x681) + x609 - x
     &611 + x614 + x615 - x616 - x617 + x625 - x627 + x632 - x635 + x637
     & + x640 + x641 + x642 + x654 + x662 - x666 + x670) - x467*x52)) + 
     &6.78155247407862d0*x228 *x341 - 8.13786296889434d0*x415*x533 - 2.0
     &7103291033722d0*x473 - 0.0638866631886015d0*x476 + x477*x479 - 0.0
     &409654285499647d0*x483 + 0.0106477771981002d0*x485 - 4.74042952861
     &565d0*x487 + x488*x497 + 1.55327468275291d0*x499*x57 + 1.553274682
     &75291d0*x500 + x502* x503 + x503*x504 - x503*x513 - 0.472898258931
     &437d0*x507 + 56.14294461557d0*x509 - 7.01786807694624d0*x510 + 35.
     &9999986138594d0*x511 - 207.78907079251d0*x512 + 246.228614519237d0
     &*x515 - 240.0d0*x516 + 7.01786807694624d0*x517 - 12.0d0*x518 - 61.
     &5571536298092d0*x520 + 96.0d0*x522 + 12.0d0* x531)
      d2E(2) = -x420 - x421 + x424 - x426 - x472 + x683*(-x171*(rho_s*x5
     &35*x591 - 2.14146774886044d0*x230*x442 - x253*x702 - 0.53536693721
     &5109d0* x427*x534 - 6.0d0*x428*x553 + x496*x545 - 6.0d0*x552*x703 
     &+ x557* x705 + x560*x704 + x587*(x211*x696 + x422*x519 + x480*x706
     & - x480 *x708 - x486*x583 - x691*x707 - 20.5190512099364d0*x695 + 
     &x697 - 4.0d0*x698 + x709) + x671*(-0.350864882184401d0*rho_s*x244*
     &x291 + 0.0799516697279215d0*rho_s*x293*x541 - x294*x465*x53 + x371
     &*( -48.0d0*rho_s*x147*x300 + rho_s*x234*x593*x595 - 24.0d0*rho_s*x
     &60 *x97*(x295 - x296) + x188*x301 - x25*x447*x731 + x311*x733*x746
     & - x333*x455 - x333*x458 - x446*x736*x763 - x447*x744 - x449*x454*
     & x736 - x451*x765 + x452*x732 - x452*x744 - x453*x454*(x757 + x758
     & + x759 - x760 - x761) + x455*x747 + x456*x457*x754 - x457*x756 - 
     &x458*x755 + x628*x736*x737 + x733*x735 + x743*(-x739 + x742) + x77
     &0) - 1.0525946465532d0*x730*x92) - x675*(0.122688256445129d0* x145
     &*x699 - x150*x729 - 0.0669208671518886d0*x427*x672 + 1.05259464655
     &32d0*x674*x723) - x682*(-0.526297323276602d0*rho_s* x379*x466 + x1
     &41*x393*(-x116*x449 + x130*x450 - x453*x773 + x456* x772 + x463 + 
     &x469 + x470*x771) + 0.789445984914903d0*x144*x62* x730 - x166*x729
     & - x414*(rho_s*x141*x677 - x109*x354*x763 + x124* x765 - x130*x732
     & - x25*x340*(-x757 - x758 - x759 + x760 + x761) - x26*x360*x754 - 
     &29.0739294853142d0*x265*x401 - x310*x735 + 1.45369647426571d0*x330
     &*x771 + 0.726848237132856d0*x331*(3.0d0* x408 - 3.0d0*x409 - x680 
     &- x681) + x333*x772 + x333*x773 + x359* x756 + x368*x731 + 21.8054
     &471139857d0*x397*x411 - x449*x633 + x596*x607 + x605*x607 - x629*x
     &737 - x630*x746 - x743*(x739 - x740 - x741) - x747*x773 + x755*x77
     &2 + x770) + 2.06952542814068d0*x532 *x723) - 1.96301210312207d0*x6
     &9*x699 - x721*(-rho_s*x568*x570* x571 - rho_s*x581*x712 + x250*x57
     &4 + 27.1408204624105d0*x250*x575 + x250*x710 - x267*x482*x711 - x4
     &32*x576*x85 - x433*x579*x718 + x435*x558*x82*x85 - x702*x713 - x70
     &2*x716 - x702*x717 - x89*( -x433*x562 + 33.420578856467d0*x433*x56
     &7 + 9.90642317524398d0* x433*x568*x572 - x433*x710 + x713*x719 + x
     &719*x720))) + 30.7785768149046d0*x181*x519 + 192.0d0*x205*x696 - 4
     &.06893148444717d0*x471*x533 - x479*x486 + 0.0204827142749824d0* x4
     &80*x686 - 17.9999993069297d0*x480*x693 - x488*x49*x691 + 0.0319433
     &315943007d0*x685 + 3.79234362289252d0*x688 - 28.071472307785d0*x69
     &2 - 246.228614519237d0*x695 + 12.0d0*x697 - 48.0d0*x698)
      d2E(3) = -rho_c*(0.0833333333333333d0*x171*(x230*x586*x780 + x253*
     &x778 + 12.0d0* x428*x703 - x49*x545 + 0.535366937215109d0*x542 - x
     &556*x705*x779 - x559*x704*x779 + x671*(-x371*(-x445 + x789) + x444
     &) + x675*( x149*x288*x780 + x443) + x682*(-x414*(-x468 - x74*x787 
     &+ x789) + x467 + 1.55214407110551d0*x532*x780) - x721*(x266*x432*x
     &435*x577 + x266*x434*x578*x718*x784 + x271*x712 + x544*x782 + x711
     &*x783 + x716*x778 + x717*x778 + x781 - x89*(x720*x778 + x781 + x78
     &2*x784 + 33.420578856467d0*x783*x784))) + 0.237021476430782d0*x182
     &*x543 - x49*x775 + x543*x777 - x583*x687 - x687*x776 + 0.854960467
     &080682d0*x694*x774)
      end subroutine


C*****************************************************************************
      pure subroutine D2ESRC_PW92_ERF_singletref_triplet(rho_c, mu, E, d
     &1E, d2E)
C*****************************************************************************
C   Implemented by E.R. Kjellgren.
C
C   Subroutine generated using Sympy 1.3
C   Generated: March 25, 2019
C*****************************************************************************
      implicit none
      real*8, intent(in) :: rho_c, mu
      real*8, intent(out) :: E, d1E(9), d2E(45)
      real*8 :: x0, x1, x2, x3, x4, x5, x6, x7, x8, x9, x10, x11, x12, x
     &13, x14, x15, x16, x17, x18, x19, x20, x21, x22, x23, x24, x25, x2
     &6, x27, x28, x29, x30, x31, x32, x33, x34, x35, x36, x37, x38, x39
     &, x40, x41, x42, x43, x44, x45, x46, x47, x48, x49, x50, x51, x52,
     & x53, x54, x55, x56, x57, x58, x59, x60, x61, x62, x63, x64, x65, 
     &x66, x67, x68, x69, x70, x71, x72, x73, x74, x75, x76, x77, x78, x
     &79, x80, x81, x82, x83, x84, x85, x86, x87, x88, x89, x90, x91, x9
     &2, x93, x94, x95, x96, x97, x98, x99, x100, x101, x102, x103, x104
     &, x105, x106, x107, x108, x109, x110, x111, x112, x113, x114, x115
     &, x116, x117, x118, x119, x120, x121, x122, x123, x124, x125, x126
     &, x127, x128, x129, x130, x131, x132, x133, x134, x135, x136, x137
     &, x138, x139, x140, x141, x142, x143, x144, x145, x146, x147, x148
     &, x149, x150, x151, x152, x153, x154, x155, x156, x157, x158, x159
     &, x160, x161, x162, x163, x164, x165, x166, x167, x168, x169, x170
     &, x171, x172, x173, x174, x175, x176, x177, x178, x179, x180, x181
     &, x182, x183, x184, x185, x186, x187, x188, x189, x190, x191, x192
     &, x193, x194, x195, x196, x197, x198, x199, x200, x201, x202, x203
     &, x204, x205, x206, x207, x208, x209, x210, x211, x212, x213, x214
     &, x215, x216, x217, x218, x219, x220, x221, x222, x223, x224, x225
     &, x226, x227, x228, x229, x230, x231, x232, x233, x234, x235, x236
     &, x237, x238, x239, x240, x241, x242, x243, x244, x245, x246, x247
     &, x248, x249, x250, x251, x252, x253, x254, x255, x256, x257
      E = 0.0d0
      d1E(:) = 0.0d0
      d2E(:) = 0.0d0
      x0 = 9.86960440108936d0
      x1 = 1d0/x0
      x2 = 2.0d0*x1
      x3 = 1.46459188756152d0
      x4 = 1d0/x3
      x5 = rho_c**(-0.333333333333333d0)*x4
      x6 = 0.194159335344114d0*x5 + 1.0d0
      x7 = 0.693147180559945d0
      x8 = -x7 + 1.0d0
      x9 = 2.14502939711103d0
      x10 = 1d0/x9
      x11 = rho_c**0.666666666666667d0
      x12 = 1d0/x11
      x13 = x10*x12
      x14 = rho_c**(-0.5d0)
      x15 = 0.564189583547756d0
      x16 = 1.41872281647967d0*x15
      x17 = 0.826307487110758d0
      x18 = rho_c**(-0.166666666666667d0)*x17
      x19 = 0.406913004517529d0*x13 + x14*x16 + 7.24010193431683d0*x18 +
     & 3.25955091942229d0*x5
      x20 = x0/x8
      x21 = x20/x19
      x22 = x8*log(0.5d0*x21 + 1.0d0)
      x23 = x22*x6
      x24 = x2*x23
      x25 = mu**2
      x26 = 0.508616435555896d0*x13*x25 + 1.0d0
      x27 = x26**(-4)
      x28 = rho_c**(-1.66666666666667d0)
      x29 = 0.0223069557172962d0*x28
      x30 = mu**5
      x31 = 1.41421356237310d0
      x32 = 0.0837224526465878d0
      x33 = exp(-0.6836107611867116d0*x5)
      x34 = 4.60115111447049d0
      x35 = 1d0/x34
      x36 = rho_c**(-1.33333333333333d0)
      x37 = x35*x36
      x38 = rho_c**1.0d0
      x39 = 1d0/x38
      x40 = 3.14159265358979d0
      x41 = 1d0/x40
      x42 = 0.009578475d0*x41
      x43 = 0.0676322201645715d0*x13 + 0.00126674656487366d0*x37 - x39*x
     &42 + 0.0188171922990732d0*x5 + 1.0d0
      x44 = x33*x43
      x45 = x32*x44
      x46 = x31*x45
      x47 = x30*x46
      x48 = x29*x47
      x49 = x25*x5
      x50 = 5.57236303607474d0*mu*x18 + 1.0d0
      x51 = 3.13317926779378d0*x49 + x50
      x52 = 1d0/x51
      x53 = mu**3
      x54 = x14*x15
      x55 = 6.76834298980504d0*x49 + x50 + 3.39260255780131d0*x53*x54
      x56 = x52*x55
      x57 = log(x56)
      x58 = x2*(x7 - 1.0d0)
      x59 = x57*x58
      x60 = rho_c**(-2.66666666666667d0)
      x61 = mu**8
      x62 = 0.00478594012208448d0*x23
      x63 = 0.133841734303777d0*x60*x61*x62
      x64 = mu**6
      x65 = x44 - 1.0d0
      x66 = 0.148394183600699d0
      x67 = x28*x66
      x68 = rho_c**(-2.0d0)
      x69 = 0.0102659822546843d0*x23
      x70 = x64*(0.0167302167879722d0*x65*x67 + 1.0525946465532d0*x68*x6
     &9)
      x71 = 0.179587122125167d0
      x72 = x44*x71
      x73 = 0.175432441092201d0*x72
      x74 = 0.825481812223657d0*x13 - 4.49737346725955d0*x5
      x75 = exp(-0.28165369188898165d0*x5)
      x76 = x40*x75
      x77 = 12.0d0*x76
      x78 = x74*x77
      x79 = 0.0524148278841779d0*x13 + 0.494402081358784d0*x5 + 1.0d0
      x80 = 1d0/x79
      x81 = -0.0259335011650457d0*x5 + 1.0d0
      x82 = x80*x81
      x83 = 2.28942848510666d0
      x84 = x34*x83
      x85 = 1.63540853354893d0*x84
      x86 = x82*x85
      x87 = x32*(x11*x86 + x38*x78)
      x88 = 0.00111534778586481d0*x28
      x89 = x31*x53
      x90 = x89*(x39*x73 + x87*x88)
      x91 = mu**4
      x92 = 0.13157433081915d0*x65
      x93 = x39*x41
      x94 = 0.022020833726518d0*x23
      x95 = 1.63540853354893d0*x84
      x96 = exp(-0.49698248213959023d0*x5)
      x97 = x96*(0.558025705063192d0*x13 - 0.352521395009435d0*x5)
      x98 = x9*x97
      x99 = 10.9027235569928d0*x98
      x100 = x66*(-x86 + x95 - x99)
      x101 = x91*(0.00139418473233101d0*rho_c**(-1.0d0)*x100 - 3.1042881
     &4221102d0*x36 *x94 - x92*x93)
      x102 = -x101 + x48 - x59 + x63 + x70 + x90
      x103 = x102*x27
      x104 = 0.0691804894611514d0*x22
      x105 = 0.129439556896076d0*x36
      x106 = x104*x105
      x107 = 2.0d0*x6
      x108 = x21 + 2.0d0
      x109 = 1d0/x108
      x110 = x19**(-2)
      x111 = x10*x28
      x112 = rho_c**(-1.5d0)
      x113 = x112*x15
      x114 = x36*x4
      x115 = rho_c**(-1.16666666666667d0)*x17
      x116 = 0.27127533634502d0*x111 + 0.709361408239834d0*x113 + 1.0865
     &1697314076d0* x114 + 1.20668365571947d0*x115
      x117 = x109*x110*x116
      x118 = x107*x117
      x119 = x26**(-5)
      x120 = x111*x119
      x121 = x101 - x48 + x59 - x63 - x70 - x90
      x122 = x121*x25
      x123 = rho_c**(-0.333333333333333d0)
      x124 = x123*x84
      x125 = 1.09027235569928d0*x124
      x126 = x12*x3
      x127 = 1.80615420514536d0*x126*x97
      x128 = 7.26848237132856d0*x123*x98
      x129 = x96*(-0.372017136708795d0*x111 + 0.117507131669812d0*x114)
      x130 = x129*x9
      x131 = 10.9027235569928d0*x11
      x132 = x130*x131
      x133 = x40*x80*x83
      x134 = 0.0141372897033723d0*x12*x133
      x135 = 1.09027235569928d0*x124*x82
      x136 = x79**(-2)
      x137 = 0.034943218589452d0*x111 + 0.164800693786261d0*x114
      x138 = x136*x137
      x139 = x138*x81
      x140 = x11*x85
      x141 = x139*x140
      x142 = x125 - x127 - x128 - x132 - x134 - x135 - x141
      x143 = 0.00139418473233101d0*x67
      x144 = 0.00232364122055169d0*x68
      x145 = x41*x68
      x146 = x106 - x118
      x147 = 1.55214407110551d0*x37
      x148 = rho_c**(-2.33333333333333d0)
      x149 = x148*x35
      x150 = -0.0225440733881905d0*x111 - 0.0031361987165122d0*x114 + 0.
     &0047892375d0* x145 - 0.000844497709915773d0*x149
      x151 = x33*(0.113935126864452d0*x114*x43 + x150)
      x152 = 0.263148661638301d0*x93
      x153 = x145*x92 + x146*x147 + 4.13905085628136d0*x148*x94 - x151*x
     &152
      x154 = 0.00168899541983155d0*x149
      x155 = x42*x68
      x156 = 0.045088146776381d0*x111
      x157 = 0.00627239743302441d0*x114
      x158 = x33*(-x154 + x155 - x156 - x157)
      x159 = x158*x71
      x160 = 0.175432441092201d0*x39
      x161 = x159*x160
      x162 = 0.122619224952946d0
      x163 = x162*x44
      x164 = 0.0399758348639607d0*x148*x163
      x165 = x68*x73
      x166 = 0.00185891297644135d0*x60*x87
      x167 = -0.550321208149104d0*x111 + 1.49912448908652d0*x114
      x168 = x38*x77
      x169 = x74*x75
      x170 = 2.14502939711103d0
      x171 = rho_c**(-0.333333333333333d0)*x170
      x172 = x134 + x135 + x141
      x173 = x32*(x167*x168 + 1.12661476755593d0*x169*x171 + x172 + x78)
     &
      x174 = x173*x88
      x175 = x25*x31
      x176 = 0.928727172679123d0*x115
      x177 = mu*x114
      x178 = 2.25611432993501d0*x177
      x179 = 1.69630127890066d0*x113*x25
      x180 = x176 + 1.04439308926459d0*x177
      x181 = -x176 - x178 - x179 + x180*x56
      x182 = 1d0/x55
      x183 = x182*x58
      x184 = x181*x183
      x185 = x158*x32
      x186 = x31*x91
      x187 = x186*x29
      x188 = x60*x66
      x189 = rho_c**(-3.0d0)
      x190 = 2.10518929310641d0*x189
      x191 = 0.0334604335759443d0*x67
      x192 = 0.101321183642338d0
      x193 = x146*x192
      x194 = 0.526297323276602d0*x68
      x195 = 0.0472353356922751d0
      x196 = mu**7
      x197 = x146*x195*x196
      x198 = 0.0669208671518886d0*x60
      x199 = rho_c**(-3.66666666666667d0)
      x200 = 0.356911291476739d0*x199
      x201 = x196*x62
      x202 = x60*x91
      x203 = 0.0571643564037363d0
      x204 = x186*x44
      x205 = x203*x204
      x206 = -0.00508309165921971d0*x189*x205 + x197*x198 + x200*x201 + 
     &0.037178259528827d0*x202*x46 + x30*(-x151*x191 + 0.027883694646620
     &3d0*x188*x65 + x190*x69 + x193*x194)
      x207 = x175*(-x161 - x164 + x165 + x166 - x174) + x184 - x185*x187
     & + x206 + x53 *(-x100*x144 + x142*x143 + x153)
      x208 = mu*x27
      x209 = 2.71262098963145d0*x120
      x210 = x182*x8
      x211 = x187*x32
      x212 = 0.172586075861435d0*x104*x148
      x213 = 0.682784063255296d0*x105*x117
      x214 = 0.129439556896076d0*x114*x117
      x215 = x10*x60
      x216 = rho_c**(-2.5d0)*x15
      x217 = x148*x4
      x218 = rho_c**(-2.16666666666667d0)*x17
      x219 = x107*x109
      x220 = x110*x219*(0.452125560575033d0*x215 + 1.06404211235975d0*x2
     &16 + 1.44868929752102d0*x217 + 1.40779759833938d0*x218)
      x221 = x107*x116**2*x20/(x108**2*x19**4)
      x222 = x116*x219*(0.542550672690039d0*x111 + x112*x16 + 2.17303394
     &628153d0*x114 + 2.41336731143894d0*x115)/x19**3
      x223 = rho_c**(-3.33333333333333d0)
      x224 = x223*x35
      x225 = x189*x41
      x226 = x33*(0.0751469112939684d0*x215 + 0.00836319657736588d0*x217
     & + 0.00394098931294027d0*x224 - 0.01915695d0*x225)
      x227 = rho_c**(-4.0d0)
      x228 = rho_c**(-4.66666666666667d0)
      x229 = x176 + x178 + x179
      x230 = mu*x180
      x231 = x230*x52
      x232 = 1.08351503479231d0*x218
      x233 = mu*x217
      x234 = x216*x25
      x235 = -x212 + x213 + x214 + x220 + x221 - x222
      x236 = 0.350864882184401d0*x189*x72
      x237 = 24.0d0*x76
      x238 = rho_c**(-1.33333333333333d0)
      x239 = x238*x84
      x240 = 0.0282745794067445d0*x12*x138*x40*x83
      x241 = 0.219734258381682d0*x217
      x242 = x140*x81
      x243 = x136*x242
      x244 = x124*x139
      x245 = x137*x242*(0.0698864371789039d0*x111 + 0.329601387572523d0*
     &x114)/x79**3
      x246 = 0.363424118566428d0*x239*x82 - x240 + x243*(0.0582386976490
     &866d0*x215 + x241) - 2.18054471139857d0*x244 - x245
      x247 = x32*x88
      x248 = x199*x65*x66
      x249 = x33*(0.455740507457808d0*x114*x150 - x189*x42 + 0.025962426
     &2672375d0* x215*x43 + 0.0375734556469842d0*x215 - 0.15191350248593
     &6d0*x217* x43 + 0.00418159828868294d0*x217 + 0.00197049465647014d0
     &*x224)
      x250 = 0.263148661638301d0*x225*x65
      x251 = x238*x98
      x252 = rho_c**(-2)
      x253 = (0.101077332976288d0*x5 + 1.0d0)*log(1.0d0 + 29.60857464321
     &6677d0/( 0.41002507094961249d0*x13 + 9.8721297225692712d0*x18 + 3.
     &2918048099450625d0*x5 + 0.76232752193529008d0*x54))
      x254 = 1.2383028969055d0*rho_c**(-2.16666666666667d0)*x17
      x255 = mu*rho_c**(-2.33333333333333d0)*x4
      x256 = 7.26848237132856d0*x239
      x257 = 0.0565491588134891d0*x133*x28 + x240 - x243*(0.058238697649
     &0866d0*x215 + x241) + 8.72217884559427d0*x244 + x245 + x256*x82
      E = -rho_c*(-x103 + x24)
      d1E(1) = -rho_c*(-x106 + x118 + 1.35631049481572d0*x120*x122 + x20
     &7*x208) + x103 - x24
      d2E(1) = -rho_c*(-2.26051749135954d0*x119*x122*x215 + 2.2994726979
     &3409d0*x121* x224*x91/x26**6 + x207*x209*x53 + x208*(mu*x181*x229*
     &x58/x55**2 - 4.52089344419912d-5*rho_c**(-4.33333333333333d0)*x204
     & - 0.0101661833184394d0*x158*x186*x189*x203 + x175*( -0.0091093036
     &3347549d0*rho_c**(-3.66666666666667d0)*x45 - 0.0799516697279215d0*
     &x148*x158*x162 + 0.350864882184401d0*x159* x68 - x160*x226*x71 + 0
     &.133252782879869d0*x163*x223 + 0.0037178259528827d0*x173*x60 - 0.0
     &0495710127051027d0*x199*x87 - x236 + x247*(-2.25322953511185d0*x16
     &7*x171*x75 - x167*x237 - x168 *(0.917202013581841d0*x215 - 1.99883
     &265211535d0*x217) - 0.751076511703951d0*x169*x170*x36 - 0.15491242
     &6780983d0*x169*x28 + x246)) + x183*(-2.0d0*x229*x231 + x230*x55*(1
     &.85745434535825d0* x115 + 2.08878617852919d0*x177)/x51**2 + x232 +
     & 3.00815243991335d0*x233 + 2.54445191835098d0*x234 - x56*(x232 + 1
     &.39252411901946d0*x233)) - x184*x231 + 0.074356519057654d0*x185* x
     &202*x31 + x195*x196*x198*x235 - x197*x200 - 0.0991420254102054d0 *
     &x199*x46*x91 - 1.30867473541471d0*x201*x228 + 0.023721094409692d0*
     &x205*x227 - x211*x226 + x30*( 0.111534778586481d0*x151*x188 - x190
     &*x193 - x191*x249 + x192*x194 *x235 - 6.31556787931922d0*x227*x69 
     &- 0.074356519057654d0*x248) + x53*(0.00619637658813783d0*x100*x189
     & - 0.00464728244110338d0*x142 *x188 + x143*(-14.5369647426571d0*x1
     &23*x130 - 3.61230841029072d0* x126*x129 - x131*x9*x96*(0.620028561
     &181324d0*x215 - 0.156676175559749d0*x217) - 0.363424118566428d0*x2
     &39 + x246 + 2.42282745710952d0*x251 - 0.299209d0*x68*x97) + 0.5262
     &97323276602d0*x145*x151 - 4.13905085628136d0*x146*x149 + x147*x235
     & - x152*x249 - 9.65778533132318d0*x223*x94 - x250)) + x212 - x213 
     &- x214 - x220 - x221 + x222) + x102*x209*x25 + 0.258879113792152d0
     &*x104*x36 - 4.0d0*x117*x6 - 2.0d0*x208*(-x175* (x161 + x164 - x165
     & - x166 + x174) - x181*x2*x210 + x206 + x211* x33*(x154 - x155 + x
     &156 + x157) + x53*(-x143*(-x125 + x127 + x128 + x132 + x172) + x14
     &4*x66*(x86 - x95 + x99) + x153))
      d2E(3) = -0.333333333333333d0*rho_c*(-0.101321996098707d0*x252*x25
     &3 + x27*(-6.0d0 *mu*x1*x210*(2.26173503853421d0*x234 + x254 + 3.00
     &815243991335d0* x255 - x56*(x254 + 1.39252411901946d0*x255)) + 4.0
     &d0*x1*x252*x57* x8 + 0.00678055584048577d0*x195*x228*x253*x61 + 0.
     &133841734303777d0*x199*x47 + 3.0d0*x64*(0.0177751651119307d0* x192
     &*x227*x253 + 0.0334604335759443d0*x248) + 3.0d0*x89*(x236 - x247*(
     &-x237*x39*x74 + x257)) + 3.0d0*x91*(-x143*( -21.8054471139857d0*x2
     &51 - x256 + x257) + 0.0524221118390615d0* x224*x253 + x250)))
      end subroutine
